EK, 0000-0002-6283-935X

The prion hypothesis states that misfolded proteins can act as infectious agents that template the misfolding and aggregation of healthy proteins to transmit a disease. Increasing evidence suggests that pathological proteins in neurodegenerative diseases adopt prion-like mechanisms and spread across the brain along anatomically connected networks. Local kinetic models of protein misfolding and global network models of protein spreading provide valuable insight into several aspects of prion-like diseases. Yet, to date, these models have not been combined to simulate how pathological proteins multiply and spread across the human brain. Here, we create an efficient and robust tool to simulate the spreading of misfolded protein using three classes of kinetic models, the Fisher-Kolmogorov model, the Heterodimer model and the Smoluchowski model. We discretize their governing equations using a human brain network model, which we represent as a weighted Laplacian graph generated from 418 brains from the Human Connectome Project. Its nodes represent the anatomic regions of interest and its edges are weighted by the mean fibre number divided by the mean fibre length between any two regions. We demonstrate that our brain network model can predict the histopathological patterns of Alzheimer's disease and capture the key characteristic features of finite-element brain models at a fraction of their computational cost: simulating the spatio-temporal evolution of aggregate size distributions across the human brain throughout a period of 40 years takes less than 7 s on a standard laptop computer. Our model has the potential to predict biomarker curves, aggregate size distributions, infection times, and the effects of therapeutic strategies including reduced production and increased clearance of misfolded protein.

A major advance in our understanding of the human brain has been the realization that the our brain is organized as a network, both at the physical and at the functional levels [1]. This quiet revolution has been made possible by the parallel development of network theory and medical imaging, in particular by the concept of small-world networks [2]. Methods originating from graph theory are now routinely used to study many aspects of brain function and the prevalent dogma is that the brain operates as an efficiently structured, modular, dynamic network with strongly connected hubs [3]. This network is optimized to rapidly transmit information, but, unfortunately, the concept of fast transport also applies to misfolded proteins that highjack the network to rapidly spread across the brain [4]. Previous work has shown that the eigenmodes of the brain network's graph Laplacian are correlated to brain atrophy in Alzheimer's disease [5], and probablistic epidemiological models have been proposed to study transference mechanisms within the network [6].

The current prevalent theory for neurodegenerative diseases is based on the prion-like paradigm [7] in which neurodegeneration is caused by the systematic invasion and conformational autocatalytic conversion of misfolded proteins [8]. In Alzheimer's disease, Amyloid beta and tau proteins are believed to act in a prion-like manner and misfold [9]. This misfolded form of the protein acts as a template on which healthy proteins misfold and grow into increasingly larger aggregates [10]. Amyloid beta is an extracellular protein that mainly spreads across the extracellular matrix [11], whereas tau is an intracellular protein that primarily propagates within the network of axonal pathways [12]. Here, we focus on tau, which spreads across the brain in a highly predictable pattern [13]: misfolded tau proteins occur first in the locus coeruleus and the transentorhinal layer from where they spread to the transentorhinal region and the proper entorhinal cortex and ultimately affect all interconnected neocortical brain regions [14]. Figure 1a illustrates the typical spatio-temporal pattern of misfolded tau protein in Alzheimer's disease inferred from histopathological observations of hundreds of human brains [15].

Understanding the progression of Alzheimer's disease is a matter of understanding the physical processes of misfolding and transport. From a modelling perspective, three approaches have been proposed to simulate the physics of neurodegeneration: (i) kinetic growth and fragmentation models to study the local interaction of aggregates of different sizes using a set of ordinary differential equations [17]; (ii) network diffusion models to study the global prion-like spreading of misfolded proteins using graph theory [5]; and (iii) reaction-diffusion-based continuum models to study the spatio-temporal evolution of patheogenic proteins using partial differential equations [18].

Figure 1b shows that continuum models with nonlinear reaction and anisotropic diffusion can accurately predict the typical pattern of tau protein misfolding in Alzheimer's disease [16]. This simulation used a Fisher-Kolmogorov model [19,20], discretized with 400 000 tetrahedral finite elements and 80 000 d.f. The continuum model displays an excellent agreement with clinical observations. However, it is computationally expensive and impractical to systematically explore a wide variety of disease and treatment scenarios. In addition, there is currently no technology to validate its predicted spreading patterns at a high enough resolution that would truly warrant a finite-element simulation with thousands of degrees of freedom. The objective of this study is therefore to create an efficient and robust simulation tool that captures the key characteristic features of pathogenic proteins in Alzheimer's disease by combining kinetic growth and fragmentation with network diffusion through a connectivity-weighted graph from the Human Connectome Project. Figure 1c suggests that-even with three orders of magnitude fewer degrees of freedom than the continuum models-our dynamic network model accurately predicts the typical spatio-temporal pattern of tau protein misfolding.

To study the kinetics of protein misfolding, we consider three popular models with different levels of complexity, the simple one-concentration Fisher-Kolmogorov model [19], the two-concentration Heterodimer model [21] and the n-concentration Smoluchowski model [22].

The simplest model to characterize protein misfolding is the Fisher-Kolmogorov model [19,20]. Initially proposed to model the spreading of a favoured gene in population dynamics, the Fisher-Kolmogorov model is now widely used to describe travelling wave solutions in ecology, physiology, combustion, crystallization, plasma physics, phase transition and biology [23]. It is based on a simple nonlinear reaction-diffusion equation for a single unknown, the misfolded protein concentration c,

where D is the diffusion tensor that characterizes global protein spreading and α characterizes the local conversion rate from the healthy to the misfolded state as illustrated in figure 2. The Fisher-Kolmogorov equation (2.1) has two steadystate solutions, an unstable steady state at c = 0 and a stable steady state at c = 1. This implies that once misfolded protein is present anywhere in the brain, c > 0, the concentration will always be repelled from the benign state, c = 0, and attracted to the misfolded state, c = 1. While the Fisher-Kolmogorov model is attractive because of its simplicity and its low computational cost, its parameter α is purely phenomenological, it provides no insight into the mechanisms of infection, and it cannot capture intermediate equilibrium states as, for example, a result of pharmocological treatment.

The simplest possible kinetic model that accounts for two configurations of the protein, the natural healthy state p and the misfolded state p, is the Heterodimer model [21]. In this model, misfolded proteins recruit healthy proteins at a rate k 11 0 , healthy proteins bind to misfolded proteins and adopt their conformation at a rate k 1 0 2 0 , and the resulting polymer fragments into infectious seeds at a rate k 2 0 2 ,

For simplicity, we collectively represent the conformational conversion from the healthy to the misfolded state as a single step through the rate constant k 12 ,

These considerations motivate a system of governing equations for the spatio-temporal evolution of the total amount of healthy and misfolded proteins p and p [24],

where D is the diffusion tensor that characterizes protein spreading, k 0 is the production rate of healthy protein, k 1 and k1 are the clearance rates of healthy and misfolded proteins, and k 12 is the conversion rate from the healthy to the misfolded state, as illustrated in figure 3.

In the initial healthy state, the healthy and misfolded protein concentrations are p 0 = k 0 /k 1 and p0 ¼ 0; in the diseased state, they converge towards p 1 ¼ k1 =k 12 and

We can simplify the Heterodimer model (2.4) by assuming that, initially, the amount of healthy protein is much larger than the amount of misfolded protein, p ) p, which implies that dp/dt ≈ 0 and r Á (D Á rp) % 0. With these assumptions, equation (2.4) provides an explicit estimate of the amount of healthy protein p,

We approximate the healthy protein concentration p using a Taylor series,

, and substitute this expression into equation (2.4),

By re-parametrizing equation (2.6) in terms of the misfolded protein concentration, c ¼ p=p max with pmax ¼ k

12 k1 =k 0 , we recover the special case of the Fisher-Kolmogorov model (2.1) for a single unknown, the misfolded protein concentration c,

Interestingly, the Fisher-Kolmogorov parameter α now takes a physical interpretation in terms of the rates of production k 0 , clearance k 1 and k1 , and conversion k 12 . While the Heterodimer model (2.4) strikes a natural balance between computational efficiency and mechanistic insight, it does not explicitly capture the size distribution of misfolded protein aggregates, their nucleation and fragmentation, and the response of the system to specific size-targeting antibodies.

To characterize the size distribution of misfolded protein aggregates, we consider the Smoluchowski model, a set of population balance equations that explicitly account for the kinetics of nucleation, aggregation, fragmentation and clearance of particles of different sizes [22]. For more than a century, the Smoluchowski model has been widely used in statistical physics to characterize processes of polymerization, coalescence of aerosols, emulsication and flocculation. It follows the n concentrations c i of particles of size i = 1, …, n and explicitly models their aggregation and fragmentation through the individual aggregation and fragmentation rates a ij and f ij with i, j = 1, …, n,

We can summarize the collective effects of aggregation and fragmentation on the concentration c i through the aggregation and fragmentation A i and F i ,

Aggregation of two smaller particles c j and c i-j creates new particles c i and removes particles c i as they aggregate with c j to larger particles c i+j . Fragmentation removes particles c i as they fragment into two smaller particles c j and c j-i and adds new particles c i from the fragmentation of larger particles c i+j into c i and c j . Taken together, the Smoluchowski model tracks the size distribution of particles through a nonlinear system of reaction-diffusion equations for the unknown concentrations c i [25],

where D i is the size-specific diffusion tensor, k 0 is the production rate, k i is the clearance rate, and A i and F i are the size-specific aggregation and clearance rates according to aggregation-fragmentation kinetics (2.9). Here, we adopt a simplification of the Smoluchowski model (2.10), the nucleated polymerization model [26,27] with a nucleus size of two and spontaneous nucleation, to model the nucleation, aggregation, and fragmentation of tau proteins in Alzheimer's disease [17] and make the following simplifying assumptions: we assume that diffusion is size-independent, D i = D; production is only possible for healthy monomers, k 01 = k 0 for i = 1, but not for misfolded particles of any other size, k 0i = 0 for i > 1; clearance occurs at k 1 for healthy monomers, is size-independent k 2 = k i for larger particles 1 < i < n; and impossible k n = 0 for the largest particle size i = n; nucleation of two monomers occurs at a nucleation rate a 11 = κ and is irreversible; aggregation of larger particles is size-independent a ij = a, but can only occur by adding single monomers for i = 1 or j = 1 and is impossible a ij = 0 otherwise; fragmentation into monomers is impossible f ij = 0 for i = 1 or j = 1, fragmentation into larger particles is size-independent f ij = f for 1 < i, j < n, and fragmentation is impossible f ij = 0 for the largest particle size i, j = n. This results in the following explicit set of equations for the concentrations of sizes 1, 2, i = 3, …, n -1, and n [28],

(2:11)

Figure 4 illustrates the Smoluchowski model for tau proteins.

In addition to the global diffusion D, this model has six local kinetic parameters, the production of healthy monomers k 0 , the clearance of healthy monomers k 1 and misfolded polymers k 2 = k i , the nucleation κ, the aggregation a and the fragmentation f. Notably, the Smoluchowski model features two distinct mechanisms to convert healthy proteins into misfolded proteins [17]: primary conversion reflected through the nucleation rate κ and secondary conversion reflected through the aggregation rate a. While the Smoluchowski model follows the size distribution of individual aggregates, allows for sizespecific transport, aggregation, fragmentation and clearance, and provides a mechanistic interpretation of the protein misfolding [29], its intrinsic disadvantage is its large number of parameters and, with it, the risk of overfitting.

A defining feature of prion-like diseases is the spreading of misfolded proteins from a small infected region along axonal fibre tracts throughout the entire brain [30]. We model this spreading as the diffusion across the brain's connectome [31], which we represent as a weighted undirected graph G with N nodes and E edges.

We extract the graph G from the tractography of diffusion tensor magnetic resonance images of 418 healthy subjects of the Human Connectome Project [32] using the Budapest

We summarize the connectivity of the graph G in terms of the degree matrix D II , a diagonal matrix that characterizes the degree of each node I, and the weighted adjacency matrix A IJ , the ratio of mean fibre number and length between nodes I and J. The difference between the degree matrix D IJ and the adjacency matrix A IJ defines the weighted graph Laplacian L IJ ,

Figure 6a illustrates the degree D II of the baseline non-weighted graph, and the degree D II of our connectivity-weighted graph G (figure 6b) along with its adjacency A IJ (figure 6c). For our connectivity-weighted graph, the degree varies between 2.1 ≤ D II ≤ 127.6, with an average degree of D II ¼ 42:8 per node, and the lowest and highest degrees in the frontal pole, shown in blue, and in the precentral gyrus, shown in red. The adjacency matrix clearly reflects the small-world architecture of our brain with strongly connected hubs within the right and left hemispheres, indicated through the lower left and upper right quadrants, and strong connections within the four lobes, indicated through the eight red regions along the diagonal. The adjacency varies between 0.01 ≤ A IJ ≤ 35.32, with an average adjacency of A IJ ¼ 1:57 per edge, and lowest and highest values between the superior parietal and the precuneus regions and between the lateral orbitofrontal and the isthmus cingulate regions. These pronounced variations in degree and adjacency confirm the general notion that the architecture of our brain resembles a smallworld network [1] in which highly connected nodes are more likely to become infected and turn into hubs of misfolded protein spreading.

We assume that the weighted Laplacian L IJ characterizes the spreading of healthy and misfolded proteins across the brain network and discretize our three kinetic models on our weighted undirected graph G. Specifically, we introduce the concentrations c I , p I , pI and c iI as global unknowns at the I = 1, …, N nodes of our graph G. This results in the discretized sets of equations for the single concentration Fisher-Kolmogorov model (2.1) with N unknowns,

for two-concentration Heterodimer model (2.4) with 2 N unknowns,

and for the n-concentration Smoluchowski model (2.10) with n × N unknowns,

We discretize our network models in time using either implicit or explicit time integration schemes to simulate the spatio-temporal evolution of misfolded proteins across the brain.

A biomarker is a global metric to characterize the evolution of neurodegeneration across the brain [36]. We calculate the biomarker abnormality as the temporal evolution of the total concentration of misfolded proteins integrated across a specific region of interest or across the brain as a whole. Biomarker abnormalities of the individual lobes provide insight into the spatio-temporal sequence of infection;

Figure 7 summarizes the biomarker abnormality in all four lobes throughout a time period of three decades as predicted by the Fisher-Kolmogorov model (3.2). For the simulation, we chose a conversion rate constant of α = 0.5 and seeded misfolded proteins by increasing the initial concentration in the entorhinal cortex to c 0 = 0.1. Our simulation uses an implicit time integration scheme with 100 times steps of Δt = 0.4 years, and runs 0.55 and 0.64 s without and with output on a standard laptop computer. Figure 1c summarizes the resulting activation sequence. We post-process the simulation to calculate the biomarker abnormality,

as the discrete sum of the misfolded protein concentration c I at the I nodes of the temporal, frontal, parietal and occipital lobes, and of the brain as a whole. All biomarker curves in figure 7 display a smooth sigmoidal form, which is in excellent agreement with brain network spreading models in general [37] and with clinical biomarker models of neurodegeneration in particular [36]. The individual biomarkers of the four lobes reveal the characteristic spreading of misfolded tau protein in Alzheimer's disease starting in the temporal lobe, shown in green, followed by the frontal, parietal, and occipital lobes, shown in red, orange, and blue. This activation sequence agrees well with the clinically observed spreading pattern [14] in figure 1a. For comparison, the dashed grey and black lines in figure 7 show the biomarker integrated over the entire brain as predicted by the continuum model [16] in figure 1b and by the Fisher-Kolmogorov network model in figure 1c. This quantitative comparison confirms that, even at a much lower spatial resolution, our network model captures the integral characteristics of continuum models for Alzheimer's disesase [38].

Figure 8 summarizes the biomarker abnormality in all four lobes, similar to figure 7, but now, instead of using the oneconcentration Fisher-Kolmogorov model (3.2), using the two-concentration Heterodimer model (3.3). For the simulation, we chose a production rate of k 0 = 1.0, clearance rates of k 1 = 0.5 and k1 ¼ 0:5, and a conversion rate of k 12 = 0.5, which results in an initial healthy concentration of p 0 = 2.0 and an initial misfolded concentration of p0 ¼ 0:0, which we increased locally in the entorhinal cortex to p0 ¼ 0:1 to seed misfolding. Our simulation uses an implicit time integration scheme with 100 times steps of Δt = 0. interpret its phenomenological rate constant α as a combination of the mechanistic rate constants k 0 , k 1 , k1 and k 12 of the Heterodimer model.

We now adopt the Heterodimer model to study the regional vulnerability of different brain regions. Specifically, we successively seed misfolded protein in all N = 83 regions, simulate the spatio-temporal spreading across the brain, calculate the resulting 83 biomarker curves, and quantify the individual infection times. Figure 9 summarizes the biomarker curves and their associated brain regions colourcoded by infection time. Misfolded proteins spread fastest when seeded in the putamen or insula with a total infection times of 20.2 years, shown in red, and slowest when seeded in the frontal pole and entorhinal region with infection times of 30.4 and 28.8 years, shown in blue. The significant variation in infection times, by more than 10 years, underlines the heterogeneity of the brain network with a few highly connected hubs [1]. These observations agree well with the hierarchical spread of epidemic outbreaks known from general network theory [39]. For comparison, the dashed grey line illustrates the lower limit of the infection time of 16.6 years, associated with a homogeneous seeding across all N = 83 regions. The mean infection time of 24.9 years on the heterogeneous network is almost exactly 50% longer. Interestingly, the entorhinal cortex, which is known as the region where misfolded tau proteins are first observed [14], is associated with the second longest in infection time. This could explain, at least in part, why tau pathology is so difficult to detect during the early stages of Alzeimer's disease [13]. The heterogeneous vulnerability of the brain network in figure 9 presents opportunities when designing treatment strategies: Reducing the local accumulation of misfolded protein in highly infectious regions such as the putamen or the insula will have a more pronounced effect on slowing down neurodegeneration than intervening in poorly connected regions such as the frontal pole. The following section focusses on different potential treatment options.

Now that we have established a solid baseline simulation for the progression of Alzheimer's disease that agrees well with clinical observations, we explore the potential of our models for simulating different treat opportunities. Two promising therapeutic strategies are currently emerging to delay or even prevent the progression of Alzheimer's disease [40]: reducing misfolding [41] and increasing clearance [42].

Even if, to date, we do not have a precise knowledge about all the model parameters, we can still perform numerical experiments to elaborate the mechanisms and time scales associated with these interventions. ) is a simple model to efficiently explore the dynamics of protein misfolding on more complex three-dimensional finite-element geometries [16] and to study the interplay of biochemical and biomechanical degeneration [38], for applications with intermediate states, we recommend using more mechanistic kinetic models like the Heterodimer model (2.4) or the Smoluchowski model (2.10). Strikingly, in the early stages of neurodegeneration, even a small reduction of misfolding can delay disease progression by several decades [41] and reduce the resting state of misfolded protein p1 ¼ k 0 = k1 À k 1 =k 12 below its untreated value, here to p1 ¼ 0:89 for k 12 = 0.45 and to p1 ¼ 0:75 for k 12 = 0.40. 14 and the α = 0.3 pattern in figure 10 show significant differences: increasing the clearance rate k1 has similar effects as decreasing the turnover rate k 12 ; it not only delays but also reduces the accumulation of misfolded protein p and with it the biomarker abnormality P. Similar to a decreased turnover, an increased clearance can delay disease progression by several decades [42] and reduce the resting state of misfolded protein p1 ¼ k 0 = k1 À k 1 =k 12 significantly below its untreated value, here to p1 ¼ 0:67 for k1 ¼ 0:60 and to p1 ¼ 0:43 for k1 ¼ 0:70. While the Heterodimer model provides valuable insight into the clearance of all misfolded proteins, it cannot predict the effect of the selective clearance of small molecules.

In this last example, we explore the interplay of nucleation, aggregation, fragmentation and the general distribution of particle size using the Smoluchowski model and highlight its advantages over the Fisher-Kolmogorov and Heterodimer models.

To illustrate the dynamics of the Smoluchowski model, we consider the biomarker spectrum of the principal moments of the size distribution, the overall concentration of misfolded proteins c i above a critical size i ≥ j,

and the overall mass of misfolded proteins i c i above a critical size i ≥ j,

where c iI denotes the concentration of proteins of length i = 1, …, n at node I = 1, …, N. Specifically, C 2 and M 2 denote the total aggregate concentration and the total aggregate mass, and c ¼ M 2 =C 2 denotes the mean aggregate length. Figure 16 illustrates the biomarker spectrum M j for a simulation with a production rate of k 0 = 1.0, clearance rates of k 1 = 0.5 for healthy monomers, k 2 = k i = 0.5 for aggregates, and k n = 0.0 for the particles of the largest size, an aggregation rate of a = 10.0, a fragmentation rate of f = 0.048, and nucleation rates of κ = 0.00016 in the entorhinal cortex and κ = 0.0 in all other regions. There is a large uncertainty about the true values for these rate constants and they may differ highly between variants [26]. Naturally, the choice of these rates will affect the sequence of events and the interplay of nucleation, aggregation, fragmentation and spread [28].

Here, we chose the parameter values such that their order of magnitude closely followed reported values in the literature [43], where the largest rate constant is the monomer production k 0 followed by the monomer clearance k 1 , the polymer clearance k i , the fragmentation f and the nucleation κ. Because of the small time constants, we now have to use a finer time discretization with 1000 times steps of Δt = 0.04 years, and we now use an explicit time integration. Our average simulation with n = 50 discrete aggregate sizes runs 1.97 and 2.20 s without and with output on a standard laptop computer; increasing the number of aggregates to n = 500 increases the simulation time to 6.96 and 16.43 seconds without and with output. Interestingly, the mass of all misfolded proteins M 2 , the dark red curve in figure 16, is almost indistinguishable from dashed black curve that highlights the biomarkers C of the Fisher-Kolmogorov model in figure 7 and P of the Heterodimer model in figure 8. However, in contrast to the Fisher-Kolmogorov and Heterodimer models, the Smoluchowski model features several competing mechanisms and time scales and allows for two distinct types of conversion, nucleation and aggregation. To better understand the dynamics of the Smoluchowski model, we explored the parameter space and learned that increasing the production k 0 , the aggregation a, or the fragmentation f accelerates and increases protein misfolding and shifts the curves in figure 16 to the left and upward; increasing the nucleation κ accelerates protein misfolding and shifts the curves to the left, but not upward; and increasing the clearance k 1 or k 2 = k i decelerates and reduces protein misfolding and shifts the curves to the right and downward.

The initial time delay between the red M 2 curve for c i ≥ 2 and the blue M 50 curve for c i = 50 illustrates the aggregation dynamics, and confirms our intuition that smaller particles have to form first to trigger the aggregation of larger particles. Intermediate aggregate sizes, 20 ≤ c i ≤ 30, highlighted through the yellow M 20 to green M 30 curves, display a small bump and increase initially, but then either partially clear or fragment into smaller aggregates.

Figure 17 illustrates the emerging aggregate size distribution and explains the dynamic features of the model: strikingly, the mean particle size, highlighted through the colourcoded dots, increases during the early stages of infection up to a mean particle size of 24.8 after 5.3 years, but then decreases gradually during the later stages towards a converged mean particle size of 15.7 after 30 years. This agrees well with the dynamics of the Smoluchowski model, which is known to predict an initial increase of the average particle size followed by a gradual decrease towards the homeostatic mean [43]. The red curve of the converged aggregate size distribution agrees well with the dashed black line of the analytical solution [43].

Figure 18 illustrates the spatio-temporal evolution of aggregates of different sizes. For illustrative purposes, rather then showing the discrete network with the individual N = 83 nodes, we colour-coded the associated 83 brain regions according to their misfolded protein concentration using the software tool Freesurfer [34]. The time lapse images show that misfolded tau proteins occur first in the locus coeruleus and transentorhinal layer from where they spread to the transentorhinal region and the proper entorhinal cortex and ultimately affect all interconnected neocortical brain regions. This spreading pattern agrees well with clinical observations [15], the continuum model of tau protein spreading [16], and our predictions with the Fisher-Kolmogorov and Heterodimer models in figure 1. Initially, the emerging aggregates are small, but they grow into progressively larger sizes as suggested by the biomarker spectrum in figure 16, and converge towards the final aggregate size distribution as indicated in figure 17. While analytical approximations exist to predict the biomarker curve and final size distribution for size-independent model parameters [43], numerical methods are necessary to predict the effect of selective parameter changes on these global readouts of the model. predicts the baseline progression of Alzheimer's disease in agreement with figure 16. Increasing the clearance of a single specific aggregate size i from k i = 0.5 to k i = 10, while keeping all other clearance rates unchanged, delays and reduces the accumulation of misfolded tau protein c i and with it the biomarker spectrum M. Clearing a specific aggregate size i at a higher rate not only affects smaller aggregate sizes through a reduced fragmentation but also larger aggregate sizes through a reduced aggregation, which, collectively, results in an overall narrower biomarker spectrum. Increasing the target size of clearance, here from c 2 shown in orange to c 6 shown in blue decelerates and reduces protein misfolding and increases narrowing of the overall spectrum. Taken together, while the Fisher-Kolmogorov model and the Heterodimer model provide valuable insight into the kinetics of protein misfolding, only the Smoluchowski model can explain the interplay between primary conversion through nucleation, secondary conversion through aggregation, and the general distribution of particle size. Although targeting a single specific size might seem rather hypothetical, we can envision therapeutic approaches that target the production or clearance of small particles below a characteristic size [44].

Despite their complexity, neurodegenerative diseases display remarkably consistent histopathological patterns. In Alzheimer's disease, these invasion patterns are highly correlated with the spreading of misfolded amyloid beta and tau proteins. Here, we modelled the spreading of tau proteins by combining misfolding kinetics and network diffusion through a connectivity-weighted graph. In our dynamic brain network model, the concentrations of healthy and misfolded protein emerge dynamically at each node and propagate across the graph through its connectivity-weighted edges. Our model correctly predicts the spatio-temporal spreading pattern of tau in Alzheimer's disease. There is currently no in vivo technology to quantify these spreading patterns longitudinally and non-invasively in humans. Our model provides a computational window into the interacting time scales and mechanisms of neurodegeneration. Its computational efficiency allows us to rapidly screen the landscape of disease-specific parameters that govern the kinetics of protein misfolding and spreading. We have demonstrated the potential of our model by simulating biomarker curves, aggregate size distributions, infection times and therapeutic intervention. A better understanding of the spreading of misfolded proteins could open new therapeutic opportunities towards blocking protein misfolding and promoting protein clearance using antibodies or small molecules. Ultimately, we envision that brain network models can help us answer some of the fundamental open questions in neurodegeneration: Why do neurodegenerative diseases progress so slowly, but yet so highly reproducibly? Can we identify early biomarkers of neurodegeneration that would allow us to interfere early? Why is neurodegeneration currently unstoppable? Can we interfere therapeutically and what would be the best time to do so? What are the roles of intra-and extracellular propagation? Can we manipulate spreading and where do we best interfere? What is the timeline of neurodegeneration? Can we predict personalized risk curves for individuals and estimate the socio-economic burden for an entire population? While we are still far from answering these questions, we believe that quantitative brain network modelling is a promising step towards identifying the key mechanisms of neurodegeneration and their roles in neurodegenerative disease.

Data accessibility. This article does not contain any additional data.

Competing interests. We declare we have no competing interests. Funding.