Sticker-and-spacer model for amyloid beta condensation and fibrillation

A major pathogenic hallmark of Alzheimer's disease is the presence of neurotoxic plaques composed of amyloid beta (Aβ) peptides in patients' brains. The pathway of plaque formation remains elusive, though some clues appear to lie in the dominant presence of Aβ1 − 42 in these plaques despite Aβ1−40 making up approximately 90% of the Aβ pool. We hypothesize that this asymmetry is driven by the hydrophobicity of the two extra amino acids that are incorporated in Aβ1−42. To investigate this hypothesis at the level of single molecules, we have developed a molecular “sticker-and-spacer lattice model” of unfolded Aβ. The model protein has a single sticker that may reversibly dimerise and elongate into semi-flexible linear chains. The growth is hampered by excluded-volume interactions that are encoded by the hydrophilic spacers but are rendered cooperative by the attractive interactions of hydrophobic spacers. For sufficiently strong hydrophobicity, the chains undergo liquid-liquid phase-separation (LLPS) into condensates that facilitate the nucleation of fibers. We find that a small fraction of Aβ1−40 in a mixture of Aβ1−40 and Aβ1−42 shifts the critical concentration for LLPS to lower values. This study provides theoretical support for the hypothesis that LLPS condensates act as a precursor for aggregation and provides an explanation for the Aβ1−42-enrichment of aggregates in terms of hydrophobic interactions.


. Introduction
Alzheimer's disease (AD), the most common cause of dementia (60-80% of cases Abeysinghe et al., 2020), and is a fatal neurodegenerative disease causing severe and devastating cognitive impairment. Age is the biggest risk factor for AD, of which, most cases are sporadic (around 95%) and occur over the age of 65 (Abeysinghe et al., Zhang et al., 2020;Zhao et al., 2020;Ayodele et al., 2021). The exact cause of AD is not fully understood and with better living conditions meaning average life expectancy in developed countries is increasing, the number of cases and the burden of neurodegenerative disease is increasing worldwide. As of 2018, an estimated 50 million people worldwide suffer from dementia, with this number expected to triple by 2050 (Patterson, 2018).
A major pathogenic hallmark in AD brains is the presence of extracellular neurotoxic plaques made up of amyloid beta (Aβ) (Stelzmann et al., 1995;Breijyeh and Karaman, 2020). Aβ is produced from amyloid precursor protein (APP) which is cleaved by α, β and γ secretases (Chen et al., 2017;Guo et al., 2020;Hampel et al., 2021). APP proteolytic cleavage can be separated into non-amyloidogenic and amyloidogenic pathways. We will focus on the amyloidogenic pathway as it is relevant to AD. First APP is cleaved by β secretase into soluble APPβ and a 99 amino acid C-terminal fragment (C99). C99 is then cleaved by γ secretase at multiple sites giving rise to Aβ of multiple lengths ranging from 39-51 amino acids (Zhang et al., 2011;Haass et al., 2012). After cleavage, Aβ is secreted from the cell where it can then form oligomers, fibrils, and finally plaques. The amyloid cascade hypothesis suggests that the deposition of Aβ into senile plaques is critical for AD pathology (Hardy and Higgins, 1992;Karran et al., 2011;Castellani et al., 2019) but the exact mechanism and change in AD brains that lead to this neurotoxic aggregation is not fully understood. Recently, this hypothesis has been called into question based on amounting evidence in disagreement. Amyloid deposition does not correlate with neuronal loss (Kametani and Hasegawa, 2018) and amyloid burden can be identified in cognitively unimpaired individuals (Arenaza-Urquijo and Vemuri, Dubois et al., 2021). Furthermore, clinical trials using therapeutics that target Aβ for degradation have been ineffective thus far (Ricciarelli and Fedele, 2017). A particular issue for anti-Aβ therapeutics is that most treatments focus on clearing insoluble aggregates (Tolar et al., 2019) with increasing evidence suggesting that pre-fibrillar soluble oligomers are orders of magnitude more toxic than fibrils and plaques (Chafekar et al., 2008;Ono et al., 2009;Sengupta et al., 2016;Huang and Liu, 2020). Unlike plaques, oligomers have been shown to impair both synaptic function and structure (Selkoe and Hardy, 2016;Kametani and Hasegawa, 2018). A recent study by Ghadami et al. (2020) has demonstrated that transthyretin is neuroprotective which is achieved by binding to Aβ oligomers, inhibiting primary and secondary nucleation without altering elongation, ephasising the role of oligomers in Aβ-mediated neurotoxicity. Consequently, understanding the role of pre-fibrillar Aβ species such as soluble oligomers and their contribution to AD pathology have gained significant interest in recent years.
Hence, in summary, while the physico-chemical properties of Aβ remains subject to debate, these properties appear to correlate with the structural properties of the aggregates. In this study, we focus our attention on the pathway of selfassembly (Michaels et al., 2020), and aim to understand how this pathway is affected by the molecular properties of Aβ, such as hydrophobicity and intrinsic disorder. We hypothesize that when Aβ 1−42 interacts with Aβ 1−40 sufficient hydrophobicity is added to enhance self-assembly, while the excluded volume interaction is limited. We will model this on a molecular level by describing the intrinsically disordered protein as chains with hydrophilic and hydrophobic segments, as well as 'sticky' building blocks that may reversibly self-ensemble into needlelike aggregates.
The aggregation of intrinsically disordered poly-peptides (IDPs) is an important topic in the biological physics of vital intracellular and extracellular processes. In particular, this class of proteins is known to undergo liquid-liquid phase separation (LLPS) to serve biological functionality, such as the (temporary) formation of membraneless organelles (Shin and Brangwynne, 2017;Jin et al., 2021), but also appears responsible to form large droplet-like precursors for the aggregation of the tau protein (Wegmann et al., 2018), which is another hallmark for AD. Large structures (possibly micellar) of approximately 50 Aβ molecules that are speculated to facilitate Aβ aggregation have also been observed (Wegmann et al., 2018). In the present study, we analyse the self-assembly of single peptides to address the current lack of physical models that explain both the formation of large agglomerates and the transition into fibrillar structures in terms of the molecular properties of the IDPs.
We argue that the general mechanism of coupled LLPS and aggregation is an example of physics that emerges from simple concepts such as (non-specific) hydrophobic interactions, excluded volume, and specific interactions. This is a hopeful scenario, as it implies that we may gain relevant molecular insight into this mechanism using strongly coarse-grained molecular models, rather than full atomistic molecular-dynamics simulations that are computationally too demanding to reach the relevant timescales. A commonly .
/fnmol. . used coarse-graining approach to capture both the polymer physics of the polypeptide and the localized formation of reversible non-covalent bonds is to develop a sticker-andspacer model. Such models were originally developed for synthetic associating polymers (Leibler et al., 1991), and were later adopted for disordered, multivalent proteins, including natural silk (Schaefer et al., 2020;Schaefer and McLeish, 2022) and scaffolding proteins in biological condensates (Shin and Brangwynne, 2017and Choi et al., 2019. Typical simulation approaches solve the Brownian dynamics of the proteins using simulation software packages such as MARTINI or LAMMPS. A major challenge, which is receiving wide attention from the modeling community (Gissinger et al., 2017;Cui et al., 2018;Raffaelli et al., 2020;Schaefer and McLeish, 2022), remains to couple the (continuum-time) conformational dynamics to model the stochastic (instantaneous) association and dissociation of reversible bonds between the molecules. To circumvent the current methodological challenges, in the present study, we will employ the fully stochastic Bond-Fluctuation Model Kremer, 1988, 1990), which by modeling the molecular dynamics on a 3D lattice successfully describes the (selfavoiding) random walk statistic of flexible and persistent chains (Bates, 2002;Feric et al., 2016;Harmon et al., 2017), LLPS, Reister et al. (2001), Reister and Müller (2003), ring-polymers (Subramanian and Shanbhag, 2009), and cross-linking reactions (Trautenberg et al., 1995). Crucially, as the dynamics are fully stochastic and discrete, they can be addressed using previously developed kinetic Monte Carlo (kMC) schemes (Lukkien et al., 1998).
In the following, we will present a simple sticker-andspacer model for the dimerisation and (linear) oligomerisation of unfolded Aβ, where the cooperativity of growth is determined both by the chain conformation (i.e., intrinsic disorder) and by the binding energies for dimerisation/nucleation, ε n , and elongation, ε e , and where hydrophobic interactions are modeled using a short-ranged interaction energy, ε H . Subsequently, we will present the bond-fluctuation model (BFM) using which we simulate the aggregation of the Aβ. In our analysis, we first investigate how the hydrophobic interactions affect the partial collapse of the molecule, and how it may lead to LLPS. We then show that these interactions facilitate dimerisation, as well as the further cooperative growth into longer oligomers. Finally, we demonstrate that increasing Aβ 1−42 relative to Aβ 1−40 promotes aggregation.
. Theory and method . . Introduction: Equilibrium statistics of self-assembly The irreversible formation of Aβ aggregates is nucleated by oligomers, which themselves are reversible and on most occasions re-dissolve into Aβ monomers (Michaels et al., 2020). The formation of oligomers can on a coarse-grained scale, be described by chemical-reaction equations of the form A 1 ⇋ A n , where n > 1 is the order of the reaction. This reaction order gives a measure for the size of the oligomers but does not provide any structural information [e.g., whether the oligomers are amorphous clusters, rings, or linear chains, and if there are parallel or anti-parallel assemblies of poly-peptide sequences (Miller et al., 2010)]. While this structural information is typically studied at the level of single amino acids and atoms, the observations on the intrinsic disorder in Aβ through its random-coil conformations through both circular dichroism (Danielsson et al., 2005) and nuclear magnetic resonance spectroscopy (Wälti et al., 2015;Roche et al., 2016), begs the question of to what extent concepts from coarse-grained polymer physics can be used to understand the difference in selfassembly behavior of Aβ 1−40 and Aβ 1−42 . The key idea that we explore is that the specific binding between substrands of the chain leads to interactions between the intrinsically disordered parts of the protein. These include both repulsive excludedvolume interactions and attractive hydrophobic interactions, ( Figure 1); their contributions grow non-linearly with the size of the oligomer.
To investigate these ideas, we will describe both species as a randomly coiled polymer that is composed of hydrophilic and hydrophobic beads, as well as one 'sticker' bead. This sticker bead forms reversible intermolecular bonds and leads to the formation of an oligomer. In general, any two oligomers of size m ≥ 1 and n ≥ 1 may reversibly assemble into an aggregate of size m + n through the chemical equilibrium equation (Martin, 1996;van der Schoot, 2005;de Greef et al., 2009).
The equilibrium constants of the reactions are non-universal and may be different for each oligomerisation step due to the formation of rings (Cates and Candau, 1990) or due to detailed phenomena, such as the parallel and anti-parallel assembly of poly-peptides. These different mechanisms may be summarised into generic schemes of (anti-)cooperative growth, where in the anti-cooperative case the nucleus is stable and only grows into larger species at high chemical potentials, whereas in the cooperative case, the nucleus is unstable and, once formed, either re-dissolves or rapidly grows into larger (stable) species. As the exact nature of the Aβ species is still under debate, we here represent these overall ideas by a simple 'dimerisation-andelongation' model (de Greef et al., 2009;Kulkarni et al., 2017), where we have a precise in silico control over the cooperativity of growth. We discuss this binding mechanism of the stickers in detail in Section 2.5. First, we will discuss the parametrisation of the protein as a randomly colloid polymer, discuss the Bond-Fluctuation model by which we simulate its equilibrium statistics, and the simple hydrophobicity model. . /fnmol. .

FIGURE
We model the growth of an oligomer through reversible bonds between bivalent stickers (red bead). The cooperativity is controlled by a di erence in the binding energy for dimerisation/nucleation (ε n ) and elongation (ε e ). The total binding energy is further a ected by repulsive excluded volume interactions between beads, as well as attractive hydrophobic interactions (ε H ) between beads within a short interaction range. The excluded volume interactions and hydrophobic interactions depend on the overlap (indicated by the dashed circles) of the intrinsically disordered sequences; the amount of overlap non-linearly increases with the length of the chains. In this study, we study how the added hydrophobic interactions and increased chain length of Aβ − compared to Aβ − a ects the self-assembly into oligomers.

FIGURE (A)
The amino acid sequence of Aβ − . Amino acids labeled red are hydrophilic, black are neutral, and blue are hydrophobic. (B) Our beads on a string model for both Aβ − and Aβ − with the N-terminus to the C-terminus going from left to right. Each bead represents amino acids. Gray beads are inert, representing the N-terminus ( -) and the turn region ( -). Blue beads are hydrophobic, representing the central hydrophobic core ( -), the second hydrophobic region ( -) and the hydrophobic C-terminus ( -/ ). Red beads are capable of dimerising with other red beads.

. . Parametrisation of unfolded Aβ
Random walk statistics of polymers typically emerge when the contour length is more than 10 times its persistence length and can be modeled using the 'Bond-Fluctuation Model' (refer to next Section). The smallest length scale of this model is the persistence length which for polypeptides is one or a few amino acids. Using these considerations, we have chosen to use 20 beads to model Aβ 1−40 and 21 beads to model Aβ 1−42 . Within our model, a bead may either be hydrophilic, hydrophobic, or may be a sticker. We have parameterised this model using the sequence of amino acids in Aβ 1−42 , (Figure 2A), which displays the hydrophilic residues in red, the neutral ones in black, and the hydrophobic ones in blue. Our model for Aβ ( Figure 2B) was based on the structure of Aβ 1−42 in Zhu et al. (2012). The first eight gray beads represent the first 16 amino acids corresponding to the hydrophilic N-terminal region. The central hydrophobic core is modeled as two blue (hydrophobic) beads representing amino acids 17-21. A red 'sticky' bead which is capable of dimerising with other red beads was added into the turn region to recapitulate the intra-and inter-chain salt bridges that can occur in this region. The remaining hydrophobic beads .
represent the C-terminal region with Aβ 1−42 containing one extra bead to account for the addition of two hydrophobic amino acids. In principle, the properties of the beads may be informed from MD simulations or structural information. However, at this high level of coarse graining the exact sequence is not expected to qualitatively alter the conclusions of our study, in which we have only observed the formation of linear oligomers (i.e., the stickers bind into supramolecular chains) and droplets, but in which no higher-order self-assembled structures such as micelles or membranes were formed (see Results). These may be expected if the hydrophobic/hydrophilic beads would have been arranged into a co-block configuration, and/or if the stickers would have been placed at the chain ends.

. . Method: Bond-fluctuation model
We model the chain conformations by placing the polymer segments on the sites of a discrete lattice, following the socalled Bond-Fluctuation Model (BFM) Kremer, 1988, 1990). This approach originates from the success of lattice models to predict phenomena such as LLPS, including in conditions near the critical point where the correlationlength diverges and mean-field models break down (Hohenberg and Halperin, 1977;Schaefer, 2018;Schaefer et al., 2019). Similarly, the lattice models have provided computationally efficient means to sample the configuration space of polymers. Early lattice models for polymers place polymer segments on a single lattice site (Binder, 1987) but underestimated the Rouse dynamics of the chains due to kinetically trapped states. This was remedied in the BFM by letting a segment occupy 8 sites on a 3D cubic lattice (4 sites on a 2D square lattice) and where the bond lengths can fluctuate from 2 to √ 10 in lattice units Kremer, 1988, 1990). This approach was proven competitive with off-lattice models both in terms of computational convenience and in physical accuracy in (and beyond) the examples mentioned in the Introduction (Trautenberg et al., 1995;Reister et al., 2001;Bates, 2002;Reister and Müller, 2003;Subramanian and Shanbhag, 2009;Choi et al., 2019).
Within the BFM, the polymer conformations and/or dynamics are modeled using kinetic Monte Carlo (kMC) time steps in which a monomer may move in 6 directions on the lattice if this does not (1) overstretch or understretch the bond between two monomers within the chain or (2) lead to double-occupied lattice sites. This leads to a list of N enabled ≤ 6N monomers possible processes that may occur during the time step, of which only one (or none) may take place. Which process, i, may take place is selected using the rate, where ν 0 is an elementary rate (typically of the order 1 − 100µs −1 ) and where E i is the change in energy. In the algorithm that we use, a process is randomly selected out of the N enabled list of processes. The process is then executed if the system energy is unchanged or decreases E i ≤ 0. However, if the energy would increase, the process is only executed with a probability p = exp(− E i /kT) and rejected otherwise; this decision is carried out using random numbers drawn using a SIMD-oriented Fast Mersenne Twister (Saito and Matsumoto, 2008). Regardless, if the process is accepted or rejected, the time is increased with t = 1/(N enabled ν 0 ). The computational efficiency of the method relies on the fact that following a Monte Carlo step during which monomer moves, only the rates of monomers in the vicinity of this monomer need to be recalculated (Lukkien et al., 1998).

. . Parametrisation: Hydrophobicity of spacers
To model the attraction between two hydrophobic monomers we follow the approach by Reister et al., and use a square-well potential (Reister et al., 2001).
where r is the distance (in units of the lattice spacing) between the two monomers, and where ε H describes non-specific (e.g., hydrophobic) interactions (Reister et al., 2001). This interaction energy enables the parametrisation of intrinsically disordered poly-peptides through the radius of gyration, and may in principle depend on conditions such the temperature and the ionic strength (Müller-Späth et al., 2010;Wuttke et al., 2014). This parametrisation is done using the Flory exponent υ, which describes the swelling of a polymer through the radius of gyration as R g = lN υ / √ 6, with l the step length between segments. In good solvent conditions, the chain is swollen due to intramolecular self-excluded volume interactions and υ = 0.588. Completely insoluble chains, described with a large value of ε H , collapse to a compact sphere with υ = 1/3. At θ conditions, the hydrophobic interactions exactly cancel the excluded-volume interactions and the chain obeys randomwalk statistics, which are characterized by υ = 1/2. To find the θ -condition, we measure the radius of gyration for chains with various chain length, N, as a function of ε H (Paul et al., 2005;Steinhauser, 2005). As R 2 g /N is independent of the chain length at the theta condition, the ε H value at which all curves intersect represents the θ condition. From Figure 3, we find that this occurs at ε H /k B T ≈ 0.27 for (homopolymer) chains with identical subunits.
As discussed in the previous section, our model for Aβ describes the protein as a copolymer with both hydrophilic .
/fnmol. . and hydrophobic units. The solid circles in Figure 3 show that the radius of gyration, R g = lN υ / √ 6 of this model polypeptide is relatively insensitive to the hydrophobic interaction parameter. To estimate the overlap concentration above which the intramolecular excluded-volume interactions are screened, we use the end-to-end distance, given by R e = lN υ ≈ 3.3 nm, where Aβ has N = 40 amino acids, and where l = 0.36 nm is the typical step length of amino acids (Müller-Späth et al., 2010;Wuttke et al., 2014). Using the molecular volume V m = (4/3)π R 3 e , we find an overlap concentration of c * = M w /(N A V m ) ≈ 100 mg/ml = 0.025 M, where we used M w = 4.5 kg/mol. In typical experiments, aggregation is observed well below the overlap concentration (e.g., Aβ 1−40 aggregates below 1 mg/ml Hortschansky et al., 2005 and Aβ 1−42 aggregates at concentrations as low as 90 nM Novo et al., 2018).
While those experimental conditions may suggest dilute conditions in which excluded-volume interactions are unimportant, inside LLPS condensates the concentration is signicantly higher, and may in fact exceed the overlap concentration. To identify LLPS in our simulations, we focus on structural coarsening phenomena such as Ostwald or Lifshitz-Slyozov-Wagner (LSW) ripening and/or Brownian coalescence (Bray, 2002). The dynamics of coarsening are typically associated with a characteristic length scale (e.g., the radius of a droplet) that grows with the one-third power of time as R * − R 0 ∝ t 1/3 with R 0 the initial length scale, which may emerge through nucleation or spinodal decomposition. In our simulations, we determine R * using the following recipe Schaefer, 2018;Schaefer et al., 2019): For a given time, t, we calculate the order-parameter field ψ(r i ), with r i the spatial coordinate of a lattice site with i the index of a lattice site. The value of ψ(r i ) is set to 1 if the site is occupied by a hydrophobic monomer, and to −1 otherwise. We then calculate the 3D Fourier transformψ(q) and obtain the structure factor S(q) = |ψ(q)| 2 , where . is the spherical/angular average. Next, we obtain the spatial correlation function C(R) as the inverse Fourier transform of S(q). As discussed previously Schaefer, 2018;Schaefer et al., 2019), numerically robust measures for the characteristic length scale R * are the first root (i.e., given by C(R * ) = 0) and the first minimum for which C(R * ) < 0 and C ′ (R * ) = 0. In our analysis, we will use both measures to assess if structural coarsening occurs in our simulations.

. . Parametrisation: Binding of stickers
As motivated in Section 2.1, we will model oligomerisation through the reversible binding of divalent stickers through a 'dimerisation-and-elongation' mechanism, where each sticker can either form a single bond or two bonds with other stickers. The dimerisation reaction is where is the equilibrium constant for nucleation, with υ the characteristic volume, H < 0 is the binding enthalpy and with S n > 0 the entropic penalty of dimerisation. In this equation, H n is controlled by the binding energy ε n > 0 and the hydrophobic interaction strength ε H . In the absence of hydrophobic interactions, H n = −ε n is exact. The entropic penalty essentially originates from an excludedvolume interaction due to a limitation of the internal degrees of freedom (DOF) of two chains that undergo dimerisation (refer to Section 3.1). In that section, we will show that these internal DOF are affected by the hydrophobicity (an increased hydrophobicity partially collapses the chain and limits the DOF prior to dimerisation), and by the concentration (above the overlap concentration the free chains have fewer DOF prior to dimerisation), such that an increasing hydrophobicity and concentration reduce the entropic penalty and enhance dimerisation. For subsequent oligomerisation steps, i.e., for n > 1, we describe the equilibrium statistics using with in which H e /k B T = −ε n +U hydrophobic +U(θ ) is composed of an oligomeration/elongation energy, ε n > 0, and hydrophobic . /fnmol. .
interactions as before. The free energy of binding is modeled as completely entropic using a square-well bending potential (BFM simulations with smoother potentials were previously discussed in the literature, refer to Weber et al., 1999;Bates, 2002;Paul et al., 2005).
for the angle, θ , between the intermolecular bonds between stickers. Here, θ = π represents fully extended bonds. Hence, for small values of θ max , the persistence length increases and the chain of stickers approaches a rigid rod. In our simulations, we set θ max = π/18 (equivalent to 10 degrees), which ensures sufficient rigidity to avoid ring formation but sufficient flexibility to avoid lattice artifacts. One of the key predictions of our simulations will be the dependence of the fraction of aggregated material on the concentration and hydrophobicity. We will interpret these findings using analytical predictions in the limit where hydrophobic interactions are absent and where the entropic penalties are constant. In this limit, there exist some known analytical predictions (van der Schoot, 2005 with σ ≡ K n /K e the so-called cooperativity factor. The concentrations of unbound Aβ is obtained from the mass balance, where ρ is the (experimentally-controlled) overall concentration of Aβ. Using the standard sum ∞ n=1 nx n = x/(1 − x) 2 for |x| < 1 and the fraction of aggregated molecules, f ≡ 1−[A 1 ]/ρ, this mass balance can be written as Martin (1996) andde Greef et al. (2009).
which provides an implicit dependence of f on ρ. This equation has three asymptotic limits of interest, namely the strongly anti-cooperative case where dimers do not grow into larger aggregates, K e → 0 (resulting in σ → ∞) leads to Martin (1996).
In the non-cooperative case, or "isodesmic" case (Smulders et al., 2010), K ≡ K e = K n (Martin, 1996), and in the strongly cooperative case, σ → 0 (van der Schoot, 2005;Smulders et al., 2010), we have f = 0 for K e ρ < 0 and for K e ρ ≥ 0. These equations show that from the anticooperative to the cooperative case, the transition from unbound to aggregated material becomes increasingly sharp. In the following, we will discuss the influence of hydrophobic and excluded-volume interactions on self-assembly in terms of the elongation constant K e and the cooperativity factor, σ .
. Results and discussion

. . Dimerisation and LLPS
To investigate how hydrophobic and steric interactions affect dimerisation, we have compared the self-assembly of chains with such interactions to the completely hydrophilic counterpart with ε H = 0, and disabled oligomerisation (i.e., ε e = 0). We have then simulated the molecular self-assembly of 100 chains in a periodic simulation box with sizes ranging from 50 × 50 × 50 to 1,000 × 1,000 × 1,000, with dimerisation energies ranging from ε n = 2 to 10k B T. The resulting fraction of dimerised chains, f , against the number density ρ (in the number of molecules per box size), is represented by the symbols in Figure 4A.
We have curve fitted each data set with fixed dimerisation energy, ε e , using the dimerisation model of Equation (12) (solid curves) with the equilibrium constant K n = K n,0 exp(ε n ). While K n,0 is explicitly independent of ε n , the curve fits yielded the apparent dependence ln K n,0 ≈ 3.0 − 1.1ε n ( Table 1). This is caused by the fact that the dimerization concentration (characterized by the inflection point of the selfassembly curve) shifts to higher concentrations (above the overlap concentration of ρ ≈ 10 −4 , see Section 2.4) for decreasing interaction strengths. Following Flory's approach to self-avoiding walk statistics, we modify the entropy for excluded-volume interactions using a mean-field description, S n = S n,0 + (∂ S n /∂ρ)ρ, and write K n = K n,00 exp(ε n + (∂[ S n /k B ]/∂ρ)ρ. Using this modification, we have been able to capture all dimerisation curves in Figure 4A with fixed K 00 = 6.85 and excluded-volume correction (∂[ S n /k B ]/∂ρ) = 1.3 · 10 3 (dashed curves). As expected, this correction is insignificant below the overlap concentration, ρ < 10 −4 , but significantly affects the self-assembly curves at higher concentrations.

FIGURE
Fraction of dimerised Aβ, f, plotted against the number density, ρ (in units of the lattice spacing). Well below the overlap concentration of ρ = − (corresponds roughly to mg/ml, refer to main text) the solution can be considered dilute, while at higher concentrations excluded volume interactions become important. (A) The dimerisation energy ε n is varied from to k B T without hydrophobic interactions (ε H = ). The solid curves are individual fits using Equation ( ) with a dimerisation constant K n = K n, exp(ε n ) with a non-constant K n, (Table ), while the dashed curves represent a simultaneous fit of all data using K n, = K n, exp[(∂[ S n /k B ]/∂ρ)ρ] with constant K = . and an entropic excluded-volume correction (∂[ S n /k B ]/∂ρ) = . · . (B) The dimerisation energy is fixed ε n = k B T, while the hydrophobicity is increased from up to ε H = . . The dashed and dotted curve are calculated using K n = . exp( + ρ) and K n = . exp( . + ρ), respectively. The data for ε H = . and . are underestimates, as ongoing LLPS enables the slow increase of aggregate sizes. The curve-fits are shown as solid lines in Figure 4. A lin-log regression yields an apparent relationship ln Kn,0 ≈ 3.0 − 1.1εn/kBT (refer to main text).
We now focus our attention on the dimerisation curve with ε n = 3k B T, which has a very low fraction of dimers below the overlap concentration but a finite fraction of dimers at higher concentrations, and increase the hydrophobicinteraction parameter, ε H , from 0 to 0.8k B T ( Figure 4B). We find that up to a value of 0.6k B T, the hydrophobicity appears to only modestly modify the equilibrium constant as K n ≈ 6.85 exp(ε n + ε H + 1300ρ); we speculate that only the central hydrophobic beads contribute to the enhancement of dimerisation. For stronger hydrophobic interactions, however, the shape of the self-assembly curve can no longer be described using the simple dimerisation model.
In fact, the datapoint for ε H ≥ 0.7 in Figure 4B is not fully converged, and the fraction of dimers increases in a slow process, as indicated in Figure 5A. This figure shows that the fraction of dimers, f , reaches a plateau at f ≈ 0.2 for times t > 100 up to t ≈ 3, 000, but then slowly increases up to f ≈ 0.4 at t = 10 5 without any sign of convergence to a higher plateau value. This process is a consequence of the slow formation of large structures, as indicated by the typical length scales in the system Figure 5B, which increase as R * − R 0 ∝ t 1/3 , which is a hallmark of LLPS, refer to Section 2.4. Here, we determined R * using the structure factor S(q, t) (top right), and its corresponding spatial correlation function, C(R, t), (bottom right). We quantified R * using the first root (C(R * , t) = 0, R 0 = 8) and the first minimum (C(R * , t) < 0 and ∂C/∂R = 0, R 0 = 12).
These findings indicate that sufficiently strong hydrophobic interactions may lead to the formation of droplets through LLPS, which inside the droplets increases the concentration of dimerising units beyond the overlap concentration (which screens the excluded-volume interactions), and provides the mass action needed to induce dimerisation.

. . Oligomerisation
Now that we have investigated the dimerisation of Aβ, we will investigate the growth of larger oligomers. For that purpose, we again first focus on the case without any hydrophobic interactions, as this enables us to isolate the impact of excludedvolume interactions on the (anti-)cooperativity of self-assembly. This also enables us to select nucleation and elongation energies Frontiers in Molecular Neuroscience frontiersin.org . /fnmol. . of interest in the subsequent simulations, in which we do switch on hydrophobic interactions. Akin to the dimerisation case, we have simulated the self-assembly of 100 Aβ 1−40 chains in box sizes ranging from 36 × 36 × 36 to 500 × 500 × 500. In these simulations, we have switched off the non-specific/hydrophobic interactions ε H = 0 and varied the nucleation energies in the range ε n = 1 − 9k B T and the elongation energies in the range ε e = 8 − 14k B T. The simulations yielded the fraction of aggregated material, f , as a function of the concentration, to which we have curve fitted the theoretical model of Equation (11). From the curvefits, we have extracted the cooperativity factor σ and the equilibrium constant K e,0 as a function of the dimerisation energy. Here, K e,0 is defined by K e = υ −1 exp(− H e + S e /k B ) ≡ K e,0 exp(ε e /k B T). (15) The results are tabulated in Table 2 and plotted in Figures 6B,C (discussed below).
A representative self-assembly curve, obtained for fixed ε e = 9k B T, is shown in Figure 6A. This panel shows the fraction of aggregated material, f , increases with an increasing number density, ρ. The sigmoidal curve becomes increasingly sharp with decreasing dimerisation energy, which, as expected (refer to Section 2.5), indicates increasing cooperativity of self-assembly. Indeed, Figure 6B shows that the logarithm of the cooperativity factor increases linearly with increasing dimerisation energy, as expected. On the other hand, the elongation constant K e,0 is in principle expected to be independent of the dimerisation energy; however, Figure 6C shows it apparently decreases with increasing dimerisation energy. We attribute this apparent dependence to the excluded-volume interactions becoming more present at higher concentrations, as we found above in the dimerisation curves of Figure 4.
To investigate the effect of non-specific hydrophobic interaction on oligomerisation, we have used the parameters ε n = 3k B T and ε e = 9k B T of a hydrophilic chain that cooperatively self-assembles at reasonably high concentrations. This ensured that low-concentration structuring due to hydrophobic interactions would not necessitate larger (and computationally expensive) box sizes. We have then used the sequences for the Aβ 1−40 and Aβ 1−42 chains and varied the hydrophobic interaction energy ε H from 0.4 to 0.7k B T. We have again simulated 100 chains in box sizes ranging from 36 × 36 × 36 to 500 × 500 × 500 lattice sites.
We have presented the results in Figures 7A,B, displaying the fraction of aggregated material f against the number density ρ for Aβ 1−40 ( Figure 7A) and Aβ 1−42 ( Figure 7B). In line with our results on dimerisation, we find that increasing hydrophobic interaction energy shifts the aggregation concentration to a lower concentration and that the transition to the aggregates state becomes sharp. For these higher hydrophobicities, we again find a slow ongoing increase of the fraction of aggregated material due to Ostwald ripening and/or rare events of fusion of condensates. After close inspection of Figures 7A,B, we observe that the transition occurs for Aβ 1−42 at lower concentrations than for Aβ 1−40 , which we attribute to a larger number of hydrophobic interactions for the longer chain. Indeed, Figure 7C shows that the mean aggregate size at ρ = 2 · 10 −4 sharply increases at ε H = 0.6 for the longer chain and at ε H = 0.65 for the shorter chain. The mean values of the aggregate size are biased by the presence of small oligomers; Figure 7D reveals the presence of aggregates with 8 or more chains.
The difference in the critical concentration for condensation of Aβ 1−40 and Aβ 1−42 , begs the question of to what extent our simple model can capture the reported observation of Aβ 1−42enhanced aggregation in mixtures of Aβ 1−40 and Aβ 1−42 (Kumar-Singh et al., 2006;Zetterberg, 2019). To investigate this, we have carried out in silico mixing experiments for two different values of hydrophobicity. First, we have chosen ε H = 0.4k B T as it previously appeared to not have led to different aggregation dynamics of Aβ 1−42 and Aβ 1−40 . Second, we have chosen ε H = 0.6k B T, as this value showed the largest difference in the aggregation propensity of the two chains (refer to Figure 7).
In Figure 8, we plot the fraction of aggregated material ( Figure 8A) and the mean aggregate size ( Figure 8B) against the fraction of Aβ 1−42 in a mixture of Aβ 1−42 and Aβ 1−40 . Having in mind that for strong hydrophobicity, structural coarsening renders the system out-of-equilibrium even at long time scales, and thermal equilibrium may never be reached (refer to Figure 5), we have plotted the results after time 4 · 10 6 ν −1 0 and after 4·10 7 ν −1 0 . As expected, for low hydrophobicity (ε H = 0.4k B T), almost no aggregation takes place and the fraction of aggregated materials is ≈ 0.1, independently of the fraction of long chains. However, for higher hydrophobicity (ε H = 0.6k B T) aggregation occurs rapidly for mixtures with 20% Aβ 1−42 (solid lines), while at long timescale mixtures with 10% Aβ 1−42 start to aggregate (dotted lines). In qualitative agreement with the experimental observations by Kuperstein et al. (2010), the critical point is strongly biased to a low fraction of Aβ 1−42 in the mixture (< 10% in our simulations). The slow dynamics are expected, as it is typical for first-order phase separation near the critical point (Hohenberg and Halperin, 1977). However, we observe another slow process beyond the critical point, namely the increase of the aggregate size at high concentrations of Aβ 1−42 (panel B). We attribute these to rare events of aggregate fusion.

. Discussion and conclusion
Our study has contributed a molecular modeling approach to address (1) the observation of Aβ 1−42 -enhanced toxic plaques associated with Alzheimer's disease, despite being the less concentrated species of Aβ, and (2) the hypothesis of condensate-precursors for fibrillation. Here, we sought a molecular explanation for the different self-assembly behavior of Aβ 1−40 and Aβ 1−42 using sticker-and-spacer models that capture the intrinsic disorder and the effect of hydrophobicity. while the two extra amino acids in Aβ 1−42 elongate the chain length, they also add hydrophobicity.
The model predicts a rich range of phenomena due to the interplay between the nucleation and elongation of aggregates with the LLPS of condensates, driven by hydrophobic or other non-specific interactions, which we have summarised in Figure 9. The high concentration inside the condensates enhances the rate of nucleation through the formation of dimers. Furthermore, at these high concentrations the excluded-volume interactions, which in dilute conditions would hamper the growth of aggregates, are screened, such that the cooperativity of oligomerisation is enhanced. Consequently, we find a sharp transition from the unbound to the aggregated state. However, the dynamics by which aggregates may form are slow close to the critical conditions for LLPS. At concentrations above the critical concentration, the growth of aggregates is slow due .
/fnmol. . to the relatively slow dynamics of Ostwald ripening and/or fusion events. A similar secondary-nucleation mechanism was discussed previously by Cohen et al. (2013) and Michaels et al. (2020), and it was found that the rate constants for primary nucleation, elongation, and secondary nucleation are 100-, 10-, and 3fold greater, respectively, for Aβ 1−42 compared to Aβ 1−40 (Meisl et al., 2014). This provides evidence that the addition of two c-terminal hydrophobic amino acids promotes aggregation propensity, in particular, the rate of primary nucleation. This is in agreement with Roche et al. (2016) who propose that the primary nucleation of Aβ is driven by non-specific hydrophobic interactions which explain the difference in the aggregation rates of Aβ 1−40 and Aβ 1−42 .
Despite its simplicity, the model also predicts that a small fraction (< 10%) of Aβ 1−42 in a mixture of Aβ 1−42 and Aβ 1−40 shifts the critical concentration for the LLPS of condensates to lower values, which in turn leads to the nucleated self-assembly of fibrillar oligomers at reduced concentrations. This finding is qualitatively consistent with Kuperstein et al. (2010) who found that even a small change from a 1:9 to a 3:7 ratio of Aβ 1−42 :Aβ 1−40 caused a dramatic change in the aggregation kinetics and toxicity of the two mixtures in vitro and in vivo. This 3:7 ratio is of particular importance as it reflects the ratio of Aβ 1−42 :Aβ 1−40 observed in familial AD patients (Scheuner et al., 1996) suggesting that maintaining a physiological ratio of Aβ 1−42 :Aβ 1−40 is of great importance and may be an effective therapeutic target.
In conclusion, we have presented a simple stickerand-spacer lattice model that captures a wide range of molecular phenomena of relevance to the literature on Aβ aggregation, which enables us to interpret those phenomena in terms of the simple concepts of nucleation-elongation .
/fnmol. .   Schematic representation of the phase diagram that we have studied in the present study. The critical concentration for oligomerisation is roughly given by ln ρ crit ∝ ε e /k B T with ε e the elongation energy. This critical concentration is a ected by excluded volume e ects near the overlap concentration and leads to an apparent dependence on the nucleation energy ε n (refer to Section . ). For a su ciently high hydrophobic interaction energy ε H , condensates form. The concentration is high within these condensates and promotes oligomer formation. The shape of this region is a qualitative estimate that acknowledges that phase separation takes place within a limited concentration range. models and hydrophobicity and polymeric excludedvolume effects at the level of coarse-grained sticky-polymer models. We hope this approach will lead to further (quantitative) refinements to the understanding of typical experimental time and length scales of Aβ aggregation in few-component in vitro and complex multi-component in vivo studies.

Data availability statement
The datasets generated/analyzed presented in this study are available on Zenodo https://doi.org/10.5281/zenodo.7053937. The simulation code is available on request.

Author contributions
CS and JC collected the data, performed the analysis, and wrote the first draft of the manuscript. CS administered the project and developed, coded the model, and the simulation algorithm. All authors conceived and designed the study, contributed to manuscript revision, read, and approved the submitted version.