Protein-protein interactions in neurodegenerative diseases: A conspiracy theory

Abstract

Neurodegenerative diseases such as Alzheimer’s or Parkinson’s are associated with the prion-like propagation and aggregation of toxic proteins. A long standing hypothesis that amyloid-beta drives Alzheimer’s disease has proven the subject of contemporary controversy; leading to new research in both the role of tau protein and its interaction with amyloid-beta. Conversely, recent work in mathematical modeling has demonstrated the relevance of nonlinear reaction-diffusion type equations to capture essential features of the disease. Such approaches have been further simplified, to network-based models, and offer researchers a powerful set of computationally tractable tools with which to investigate neurodegenerative disease dynamics. Here, we propose a novel, coupled network-based model for a two-protein system that includes an enzymatic interaction term alongside a simple model of aggregate transneuronal damage. We apply this theoretical model to test the possible interactions between tau proteins and amyloid-beta and study the resulting coupled behavior between toxic protein clearance and proteopathic phenomenology. Our analysis reveals ways in which amyloid-beta and tau proteins may conspire with each other to enhance the nucleation and propagation of different diseases, thus shedding new light on the importance of protein clearance and protein interaction mechanisms in prion-like models of neurodegenerative disease.

About the Authors:

Travis B. Thompson

Contributed equally to this work with: Travis B. Thompson, Pavanjit Chaggar, Ellen Kuhl, Alain Goriely

Roles Conceptualization, Formal analysis, Investigation, Methodology, Software, Visualization, Writing – original draft, Writing – review & editing

* E-mail: [email protected]

Affiliation: Mathematical Institute, University of Oxford, Oxford, United Kingdom

ORCID logo https://orcid.org/0000-0001-8696-4115

Pavanjit Chaggar

Contributed equally to this work with: Travis B. Thompson, Pavanjit Chaggar, Ellen Kuhl, Alain Goriely

Roles Data curation, Investigation, Software, Visualization, Writing – original draft, Writing – review & editing

Affiliation: Mathematical Institute, University of Oxford, Oxford, United Kingdom

ORCID logo https://orcid.org/0000-0002-8395-3047

Ellen Kuhl

Contributed equally to this work with: Travis B. Thompson, Pavanjit Chaggar, Ellen Kuhl, Alain Goriely

Roles Conceptualization, Funding acquisition, Supervision, Writing – review & editing

¶‡ Senior Author.

Affiliation: Living Matter Laboratory, Stanford University, Stanford, California, USA

ORCID logo https://orcid.org/0000-0002-6283-935X

Alain Goriely

Contributed equally to this work with: Travis B. Thompson, Pavanjit Chaggar, Ellen Kuhl, Alain Goriely

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Software, Supervision, Visualization, Writing – original draft, Writing – review & editing

¶‡ Senior Author.

Affiliation: Mathematical Institute, University of Oxford, Oxford, United Kingdom

ORCID logo https://orcid.org/0000-0002-6436-8483

for the Alzheimer’s Disease Neuroimaging Initiative

¶Data used in preparation of this article were obtained from the Alzheimer’s Disease Neuroimaging Initiative (ADNI) database (adni.loni.usc.edu). As such, the investigators within the ADNI contributed to the design and implementation of ADNI and/or provided data but did not participate in analysis or writing of this report.

Introduction

Neurodegenerative diseases such as Alzheimer’s (AD) or Parkinson’s (AD) are associated with the propagation and aggregation of toxic proteins. In the case of AD, it was Alzheimer himself who showed the importance of both amyloid-β (Aβ) plaques and tau-protein (τP) neurofibrillary tangles (NFT) in what he called the “disease of forgetfulness” [1, 2]. These two proteins are very different. Aβ forms extracellular aggregates and plaques whereas τP are intracellular proteins involved in the stabilization of axons by cross-linking microtubules that can form large disorganized tangles [3, 4]. Since the early 90’s, when it was first formulated, the “amyloid cascade hypothesis” has dominated the search for cures and treatments [5, 6]. According to this hypothesis, an imbalance between production and clearance of Aβ42 and other Aβ peptides is not only an early indicator of the disease but the causing factor for its initiation, progression, and pathogenesis [7]. However, the repeated failures of large clinical trials focussing on the reduction of Aβ plaques has led many researchers to question the amyloid hypothesis and argue for the possible importance of other mechanisms.

One obvious alternative is that τP plays a more prominent role than the amyloid hypothesis suggests. The τP are usually considered as secondary agents in the disease despite the fact that (1) other τP-related diseases (tauopathies), such as frontotemporal lobar degeneration, are mostly dominated by τP spreading [8]; (2) brain atrophy in AD is directly correlated with large concentrations of NFT [9, 10]; (3) τP distribution determines disease staging [11]; (4) lowering τP levels prevent neuronal loss [12]; (5) τP reduces neural activity and is the main factor associated with cognitive decline [13]. These findings may explain the relative lack of clinical improvements after Aβ suppression and the debate between the relative importance of Aβ proteopathy and τP tauopathy in AD [14]. Furthermore, the similarity in mechanism and progression between prion diseases [15] and classical neurodegenerative diseases led to the formulation of the “prion-like hypothesis” [16, 17, 18, 19, 20, 21] stating that all these protein-related degenerative diseases are characterized by the progressive spreading and autocatalytic amplification of abnormal proteinaceous assemblies through axonal pathways [22].

Since so many cellular mechanisms are poorly understood in vivo, the relative importance of different groups of toxic proteins and their possible interactions have not been established. In particular, both τP and Aβ depend upon and modify the cellular environment [16]. Yet, in recent years a number of studies have linked these two anomalous proteins [23] and raised the possibility that protein-protein interactions in neurodegenerative diseases are a key to understanding both spreading and toxicity [24, 25]. According to Walker, for AD “the amyloid-β-τ nexus is central to disease-specific diagnosis, prevention and treatment” [14].

Specifically, the following crucial observations have been made in AD: (1) tangles in the cortex rarely occur without Aβ plaques [12]; (2) the presence of Aβ plaques accelerates both the formation of τP aggregates [26] and the interneuronal transfer [27] of τP; (3) the presence of τP induces blood vessel abnormalities [28] and induces neuroinflammation through micro- and astro-glial activation [29]; (4) the presence of Aβ can induce the hyperphosphorylation of τP and the creation of toxic τP seeds by disturbing cell signaling via oxidative stress or through plaque-associated cells (such as microglia) or molecules [30, 31, 32, 33]; (5) Aβ and toxic τP target different cellular components, and doing so amplify each other’s toxic effects [23]; (6) τP mediates Aβ toxicity as a reduction of τP prevents Aβ-induced defects in axonal transport [23]; (7) perhaps more anecdotal, it has also been argued that the lack of clear evidence of dementia in non-human primates, despite the presence of Aβ plaques, could be due to a difference in Aβ-τ interactions in these species [34].

From these observation, we extract three crucial modes of interaction:

1. M1: The seeding of new toxic τP is enhanced by the presence of Aβ.

2. M2: The toxicity of Aβ depends on the presence of τP.

3. M3: Aβ and τP enhance each other’s toxicity.

Here, our goal is twofold: first to develop modeling and computational tools to study protein-protein interactions at the brain-organ level and second to test the relative effect of these interactions by direct simulation. Typical approaches for organ-size simulation of dementia progression [35] take the form of either continuous models formulated in terms of anisotropic reaction-diffusion equations [36, 37, 38], or discrete systems on the brain’s connectome network. The discrete approach can be further divided into pure-diffusion linear models [39, 40, 41, 42, 43], probabilistic models [44, 45, 46], or deterministic models [47, 48].

A primary result, of interest to the computational biology community, for the current work will be to show: that non-trivial interactions between Aβ and τP can be realized with relatively simple deterministic models and couplings; and that these interactions can lead to effects with physiological interpretations in neurological disease modeling. Moreover, the mathematical analysis will highlight that clearance mechanisms play a key role in destabilizing the system towards proteopathy. We will therefore select the simplest possible, deterministic, protein kinetic model, including a bulk clearance term, that allows for the expression of both a healthy and toxic regime for a single protein; the heterodimer model [49, 50]. One such system will be defined for Aβ, one for τP and the two heterodimer systems are coupled with a single balanced interaction term. We augment the model by adding an stand-alone, first order equation for damage evolution; this equation expresses the deleterious effects of Aβ, τP and their interactions. Our general approach, following [48] is to study some of the key properties of this continuous model before discretizing it on a network and solving it numerically on the brain’s connectome graph.

Results

Network model dynamics

We have established the properties of our system of equations in the homogeneous case and in one-dimension (Methods, A network mathematical model). The study has lead to the identification of two fundamental disease propagation modes depending on the parameters: the primary tauopathy where toxic τP states can exist independently from the Aβ concentration, but are enhanced by its presence; and the secondary tauopathy where the presence of toxic τP is slaved to the existence of toxic Aβ. We can use this analysis as a guide to the simulation of the full network equation. Eqs (8)–(12) were discretized on the reference connectome [59] (Methods, A network mathematical model) using CVODE as part of the SUNDIALS nonlinear ODE solver library [63] in addition to KLU [64] as part of the SuiteSparse [65] linear algebra library. Snapshots of the dynamics are shown in subsequent figures, but full videos can be found in the supplementary material.

As a way to systematically test the validity of our computational platform, we have performed two main tests. First, we reproduced the homogeneous states in the full network and second, we reproduce the transition between homogeneous states. Both tests are detailed in S1 Appendix in addition to a discussion regarding a choice of hypothetical, non-clinical parameters for illustration purposes; S1 Appendix contains a full discussion on the specific numeric values of the stationary states corresponding to the choices of of Table 1.

[Figure omitted. See PDF.]

Table 1. Primary tauopathy model parameters.

https://doi.org/10.1371/journal.pcbi.1008267.t001

Front dynamics on networks.

Propagating front solutions for the system of partial differential Eq (3) were considered, via linearization around the healthy state and reduction to one spatial dimension, (c.f. Methods, Front propagation), Propagating fronts represent fundamental modes of disease pathology dynamics that can also be realized by the network model of (8)–(12) as we now demonstrate. We consider two different network for front propagation. First, a three-dimensional regular cubic lattice with nx = 30 nodes in the x-direction ny = 6 nodes in the y-direction and nz = 3 nodes in the z-direction, spaced equally at unit length. Second, we use the physiological brain connectome domain discussed in the Methods section (A network mathematical model), but we choose initial conditions on two sides of the brain to illustrate the front dynamics. In the next section we will consider the same domain but with realistic initial conditions.

The first example is that of primary tauopathy corresponding to the parameters of Table 1. Formulas for the steady states, for primary tauopathy, are listed in the Methods (An Analysis of the continuous model, Stationary points) section.

Primary tauopathy.

Primary tauopathy, synthetic domain.

We set all nodes to the healthy state and perturb the initial condition of the left-hand nodes 0 ≤ x ≤ 4 by adding a 5% concentration () of toxic Aβ. We perturb the initial condition of the right-hand nodes 25 ≤ x ≤ 29 by adding a 5% concentration () of toxic τP. As expected, we see the toxic Aβ concentration achieve the theoretical maximum, permitted by the parameters, of while toxic τP first achieves the maximum associated with and, upon mixing with Aβ, achieves the fully toxic state value . The color scale of Fig 1 was chosen to accentuate the interaction.

[Figure omitted. See PDF.]

Fig 1. Front propagation in primary tauopathy; synthetic rectangular domain.

Each subfigure consists of a toxic Aβ concentration distribution (top left), toxic τP concentration distribution (bottom left) and a plot (solid line: Aβ, dashed line: τP) of the concentration level on the x−axis. Dark blue indicates the minimum concentration of c = 0.0 while bright red indicates the maximum of c = 0.5. See the Methods section (Front propagation) for a comparison to theory. (See also: supplementary S1 Video).

https://doi.org/10.1371/journal.pcbi.1008267.g001

Primary tauopathy, brain connectome.

Simulation of disease front propagation was then carried out using the physiological connectome (Methods, A network mathematical model). The seeding sites selected for toxic Aβ and toxic τP are the right supramarginal gyrus and left supramarginal gyrus respectively; these seeding sites provide a direct analogy, when the brain connectome is viewed from the frontal lobe, with Fig 1. Fig 2 depicts time instances qualitatively reflecting, in one-to-one correspondence, the stages of the synthetic domain computation of Fig 1. A horizontal slice, at the plane of the supramarginal gyri, of the brain connectome is used to maximally expose the front propagation dynamics. The impact of brain connectome cross-connectivity is evident in the stages depicted in Fig 2. In particular, when the Aβ and τP wavefronts first meet they do so in several locations. This is due to the left-right hemispheric connectivity; both direct nodal connectivity and vis-a-vis propagation in the coronal plane.

[Figure omitted. See PDF.]

Fig 2. Front propagation in primary tauopathy; brain connectome.

Each subfigure consists of a toxic Aβ concentration distribution (subfigure left) besides a toxic τP concentration distribution (subfigure right). Dark blue indicates the minimum concentration of c = 0.0 while bright red indicates the maximum of c = 0.5. (See also: supplementary S2 Video and supplementary fle S2 Data).

https://doi.org/10.1371/journal.pcbi.1008267.g002

Secondary tauopathy.

The parameters for the at-risk secondary tauopathy patient are those of Table 1 with two exceptions; first, as usual for secondary tauopathy, we take b2 = 0.75 and second we take b3 = 3.0. We have increased b3 to facilitate the comparison with front propagation theory, for secondary tauopathy, discussed in the Methods section.

Secondary tauopathy, synthetic domain.

Secondary tauopathy consists of all stationary states except for the toxic τP–healthy Aβ state; i.e. is not included. The stationary point depends on b3; with the parameters above we have(1)while the other two secondary tauopathy stationary points, c.f. (14) and (15), coincide with their values for primary tauopathy. The initial value at all nodes are first set to the healthy state. A 5% perturbation in concentration is then added to the toxic Aβ initial value for the nodes 0 ≤ x ≤ 4 and a perturbation of 1 × 10−9%, i.e. 1 × 10−11, is added to the toxic τP initial value for the nodes 0 ≤ x ≤ 14.

As expected: the initial toxic Aβ wavefront achieves its theoretical maximum of ; c.f. Fig 3 and the front propagation discussion, for secondary taopathy, in the Methods section. The toxic τP wave takes on detectable concentration levels at the point when the Aβ wave reaches the halfway mark in the rectangular domain. The toxic τP state connects, immediately, to the theoretical maximum of the toxic τP–toxic Aβ stationary state value of and quickly proceeds to catch up to the Aβ wavefront.

[Figure omitted. See PDF.]

Fig 3. Front propagation in secondary tauopathy; rectangular domain.

Each subfigure consists of a toxic Aβ concentration distribution (top left), toxic τP concentration distribution (bottom left) and a plot (solid line: Aβ, dashed line: τP) of the concentration level on the x−axis. Dark blue indicates the minimum concentration of c = 0.0 while bright red indicates the maximum of c = 0.5 for toxic Aβ and for toxic τP. See the Methods section (Front propagation) for a comparison to theory. (See also: supplementary S3 Video).

https://doi.org/10.1371/journal.pcbi.1008267.g003

We tested the time of appearance and saturation of the toxic τP wave front as a function of the interaction parameter b3. Plots for four values of b3 are shown in Fig 4 where the y-axis signifies the maximal toxic τP concentration obtained, over all nodes, with respect to the maximum concentration for that value of b3 (c.f. (1)). Fig 4 highlights the important, and patient-specific, role that b3 may play in further efforts to deploy (8)–(11) for the modeling of Alzheimer’s disease. In particular values of b3 ≈ 1 do lead to the development of tauopathy; however, this development emerges significantly later than for higher values of this interaction parameter. Clinically, such a value of b3 could correspond to a patient who, at the time of death, presents significant amyloid plaques but negligible, or undetectable, levels of neurofibrillary tau tangles.

[Figure omitted. See PDF.]

Fig 4. The onset effect due to b3 in secondary tauopathy.

Saturation % (y-axis) vs Simulation time (x-axis).

https://doi.org/10.1371/journal.pcbi.1008267.g004

Secondary tauopathy, brain connectome.

We also simulated secondary tauopathy dynamics on the physiological brain connectome (Methods, A network mathematical model). A 5% toxic Aβ perturbation from the healthy state was seeded at the site of the left supramarginal gyrus; all nodes of the left hemisphere were then seeded with an additional 1 × 10−9% concentration of toxic τP. Snapshots of the evolution is shown in Fig 5. As indicated above we have b3 = 3 for comparison with Fig 3 and the theory of propagating fronts as illustrated in the Methods section. A detail of particular interest is that, even though the entire left hemisphere was seeded uniformly with toxic τP, the toxic τP wave follows the same anisotropic infection pathway, from the left supramarginal gyrus, as the toxic Aβ front propagation. This implies that latent development of tauopathy, in this regime, is heavily influenced by Aβ pathology history.

[Figure omitted. See PDF.]

Fig 5. Front propagation in secondary tauopathy; brain connectome.

Each subfigure consists of a toxic Aβ concentration distribution (subfigure left) besides a toxic τP concentration distribution (subfigure right). Dark blue indicates the minimum concentration of c = 0.0 while bright red indicates the maximum of c = 0.5 for toxic Aβ and for toxic τP. (See also: supplementary S4 Video and supplementary S2 Data).

https://doi.org/10.1371/journal.pcbi.1008267.g005

Application to neurodegenerative disease modeling

We have shown in the previous section that the overall phenomenology obtained from the dynamic evolution of the continuous model in one-dimension (Methods, A network mathematical model, Front Propagation) is recovered within the discrete network setting. We can therefore use the network model and our primary classification to study the interaction of proteins in the brain. Here, we apply (8)–(11) to a computational case inspired by Alzheimer’s disease. In particular we consider seeding sites, for toxic Aβ and toxic τP, commensurate with [11, 66, 67, 47] Alzheimer’s disease staging. Alzheimer’s disease is a complex multiscale phenomena; a uniform parameter regime, throughout all brain regions, is unlikely to accurately reflect a patient’s real disease progression. Nevertheless, for this early investigation, we will first consider the simple uniform parameters, of the model’s primary and secondary tauopathy regimes, as discussed in the section on network model dynamics. In addition we briefly consider the evolution of the coupled neuronal damage term, given by (12), and the effect of the coefficients therein. We shall also select the diffusion constants, ρ of (7), to be unity for (8)–(11).

A simplified model of Alzheimer’s disease proteopathy.

Alzheimer’s associated amyloid deposition begins [18, 47, 66, 67] in the temporobasal and frontomedial regions. Tau staging, in Alzheimer’s disease, follows the Braak tau pathway [11] and begins in the locus coeruleus and transentorhinal layer [18, 47, 67]. These seeding sites, used throughout this section, are shown in Fig 6. The temporobasal and frontomedial regions for toxic Aβ seeding are highlighted in red on the left while the locus coeruleus (in the brain stem) and transentorhinal associated regions, for toxic τP staging, are highlighted red on the right.

[Figure omitted. See PDF.]

Fig 6. Simulated seeding sites for a model of Alzheimer’s disease.

Toxic Aβ (left) and toxic τP (right) seeding sites.

https://doi.org/10.1371/journal.pcbi.1008267.g006

The regimes of primary and secondary tauopathy will first be considered on the whole brain connectome with globally-constant synthetic parameters. We will observe several characteristic traits of these modalities and also note the similarity between these pure states and to a qualitative three-stage progression [18] of protein lesions, typical of Alzheimer’s disease, as inferred from post-mortem analyses; the progression pattern is illustrated in Fig 7. In the subsection titled ‘a mixed model comparison to Alzheimer’s diseased patient data’ we consider the case of mixed regional modalities; i.e. a mixture of primary and secondary tauopathy connectome regions. We illustrate that the model can manifest canonical features of positron emission tomography (PET) SUVR intensities characteristic of Alzheimer’s disease (c.f. for instance [68, 69]). In particular: we will compare the results of a mixed-mode simulation with a cross sectional Alzheimer’s patient cohort dataset procured from the Alzheimer’s Disease Neuroimaging Initiative (ADNI) database.

[Figure omitted. See PDF.]

Fig 7. Characteristic progression of of Aβ and τP lesions.

3-stage Aβ (top) progression and τP NFT (bottom) progression.

https://doi.org/10.1371/journal.pcbi.1008267.g007

Alzheimer’s disease, primary tauopathy.

All nodes in the connectome were first set to the healthy, but susceptible, primary tauopathy patient stateThe temporobasal and frontomedial Aβ seeding sites, consisting of fifty-three nodes, were each seeded with a toxic amyloid concentration of 0.189%; thus the brain-wide toxic Aβ concentration represents a 1% concentration deviation from healthy. Similarly, the locus coeruleus and transentorhinal nodes were seeded with an aggregate perturbation of 1% toxic τP.

Fig 8a shows the average brain-wide concentration for all four protein populations for the primary tauopathy patient (c.f. Table 1) with interaction term b3 = 1. As we observed previously, in (1), the value of b3 directly informs the saturation τP concentration, of , for the disease. Fig 8b shows the evolution of the toxic τP burden for various b3. For each value of b3 the toxic τP invasion window was computed as the difference in time between the appearance of a global 1% toxic τP concentration to the simulation time where the maximum was reached. We performed a least squares fit and found that the invasion window, for primary tauopathy, decreases exponentially with an increase in coupling strength (b3) between toxic Aβ and toxic τP. Fig 8c shows the result.

[Figure omitted. See PDF.]

Fig 8. Protein-protein interaction in primary tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.g008

This result suggests that the dynamics of toxic protein evolution is highly sensitive to the coupling between Aβ and τP: Toxic Aβ accelerates, in a nonlinear fashion, the way toxic τP emerges across the brain. Acceleration of toxic τP progression due to the presence of toxic Aβ has also been observed in mouse models of Alzheimer’s disease [27]. Consulting longitudinal tau PET studies, in combination with amyloid-beta data from a public database, could provide an estimation of b3 in the primary tauopathy model.

The toxic load progression of the susceptible primary tauopathy patient is shown in Fig 9 at five equidistant time points throughout the invasion window. To facilitate a comparison with Fig 7: a sagittal view of the progression, of each toxic agent, is presented; directly below is an opacity-exaggerated view wherein regional opacity is proportional to the agent’s regional toxic load. Comparison with Fig 7 suggests reasonable qualitative agreement; thus warranting further study of physically relevant parameters with a view towards real clinical applications.

[Figure omitted. See PDF.]

Fig 9. Toxic proteopathy progression dynamics in the primary tauopathy patient.

Toxic Aβ (top row) and opacity exaggerated toxic Aβ progression (second row); Toxic τP (third row) and opacity exaggerated toxic τP progression (last row). Color scale is identical to Fig 1. (See also: supplementary S5 Video, supplementary S6 Video and supplementary S2 Data).

https://doi.org/10.1371/journal.pcbi.1008267.g009

Alzheimer’s disease, secondary tauopathy.

All nodes were set to the healthy, but susceptible, patient state corresponding to the susceptible secondary tauopathy patient parameters (Table 1 with b2 = 0.75). In addition, for a baseline secondary tauopathy case, we follow the secondary tauopathy approach discussed previously and select the interaction parameter of b3 = 3.0; the fully invaded secondary tauopathy state values are therefore (1). Seeding patterns for both Aβ and τP are identical to the case of primary tauopathy discussed above.

Fig 10a shows the average brain-wide concentration for all four protein populations of the secondary tauopathy patient with baseline interaction term b3 = 3. As in the case of primary tauopathy we investigate the effect of b3 on toxic load and invasion window by considering a value range four times smaller to four times larger than the baseline b3 = 3 case. Toxic load curves are shown in Fig 10b while invasion windows are shown in Fig 10c.

[Figure omitted. See PDF.]

Fig 10. Protein-protein interaction in secondary tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.g010

Interestingly, we see distinct differences in comparison with the primary tauopathy case (c.f. Fig 8a–8c). More specifically, in primary tauopathy it is evident (Fig 8b) that the disease onset is only slightly affected by varying the interaction parameter b3; for secondary tauopathy, in contrast, b3 has a profound effect on disease onset latency. Moreover, the invasion window variation with b3 for secondary tauopathy is more complex than that of primary tauopathy. Fig 10c shows that the invasion window duration initially decreases exponentially with b3 but then appears to increase logarithmically for b3 ≥ 3. Analyzing the invasion window start time and end time separately shows a clear, but separate, exponential decay pattern versus b3. Fig 11 shows the least-squares exponential fit to the invasion start and end times.

[Figure omitted. See PDF.]

Fig 11. Prodromal window variations with b3, secondary tauopathy.

Invasion starting (left) and ending (right) time vs. b3.

https://doi.org/10.1371/journal.pcbi.1008267.g011

As in the primary tauopathy case we now consider characteristic toxic load progression for secondary tauopathy. The Aβ progression is identical to that shown in Fig 9 (top two rows). This is expected as only the τP portion of the system has been modified with respect to the primary tauopathy regime (S1 Appendix). The τP secondary tauopathy progression is shown, in Fig 12, at equally spaced simulation times through the invasion window. Qualitatively, the progression of secondary tauopathy also reflects the characteristic post-mortem progression of Fig 7.

[Figure omitted. See PDF.]

Fig 12. Toxic τP progression dynamics in the secondary tauopathy patient.

Toxic τP (first row) and opacity exaggerated toxic τP progression (second row). Color scale is identical to the τP case of Fig 3. (See also: supplementary S5 Video, supplementary S7 Video and supplementary S2 Data).

https://doi.org/10.1371/journal.pcbi.1008267.g012

A mixed model comparison to Alzheimer’s diseased patient data.

Thus far, we have considered, respectively, the general features of the modalities of primary versus secondary tauopathy; illustrated with synthetic, globally constant parameters. We have observed several interesting facets of these two disease states. For instance, regions in a state of primary tauopathy can develop Aβ and τP proteopathy separately; the Aβ interaction parameter, b3, does not alter the onset of tauopathy but does modulate the regional concentration. Conversely, in secondary tauopathy the presence of Aβ pathology is necessary for τP pathology and the interaction parameter, b3, modulates both the latency and the intensity of the regional pathology. We have also seen that proteopathy progression in pure models, e.g. where all regions have the same primary or secondary tauopathy parameters, bears a notable resemblance (Fig 7) to post-mortem progression of protein lesions [18]. However, PET imaging studies of Aβ and τP radiotracer uptake tell a more nuanced story. For instance, in Alzheimer’s disease the distribution of ([18F]flortaucipir and [18F]THK-5117, among others) PET-τP SUVR intensities are distinctly biased [68, 69] towards the temporal and parietal regions of the brain; a feature that we do not see in Fig 9 or Fig 12.

In order to demonstrate that the model of (8)–(11) can reproduce salient features of Aβ and τP SUVR uptake in patients diagnosed with Alzheimer’s disease: we now compare a mixed-modality simulation with a cross-sectional study of Alzheimer’s disease patient data. Sample data for model comparison was procured from the the Alzheimer’s Disease Neuroimaging Initiative (ADNI) database. We first queried the ADNI database to locate AD-diagnosed subjects between the age of 70 and 90 who had at least one τP PET (18F-AV1451, flortaucipir) scan. The returned results consisted of 41 patients. These initial patient candidate IDs were then checked for a structural T1 weighted MRI within a maximum of one year of a tau PET scan; patients without any sMRI, those with a poor quality sMRI, or those without an sMRI within a year of the tau PET scan were discarded. The resulting cohort consisted of 38 patients (25 male and 13 female) with mean age 78.4, a standard deviation of 5.14 years, and a male-to-female ratio of 1.92. Patients who met the age, τP PET scan, and acceptable quality sMRI within one year of the PET scan were not abundant in the ADNI database. Due to this the τP group was selected, first, to maximise the number of candidates in the group.

The ADNI database contains quite a generous number of patients with Aβ (18F-AV45, florbetapir) PET scans. We next queried the database to locate AD-diagnosed subjects between 70 and 90 who had an Aβ (18F-AV45) PET scan in addition to a structural T1 weighted MRI within one year. The result of this search was in excess of 100 unique patient IDs; from these results we selected an initial candidate group of 82 unique patients IDs. The 82-candidate group was further pruned to create a list of 48 subjects whose age and sex characteristics closely resembled that of the τP PET group. Finally, the 48 candidate Aβ PET group was narrowed down: first, subjects with an unacceptable or low-resolution sMRI were removed. We then removed the minimum number of candidates required to provide as close a match as possible to the mean age, standard deviation and male-to-female ratio of the τP PET group. The resulting AD cohort for AB consisted of 42 patients with mean age 78.4, a standard deviation of 5.1 years, and a male-to-female ratio of 2.0. The two groups are succinctly summarised in Table 2.

[Figure omitted. See PDF.]

Table 2. Adni patient group statistics.

https://doi.org/10.1371/journal.pcbi.1008267.t002

Patient data was then processed through a semi-automated, scripted software pipeline for general connectome-graph based imaging and analysis of clinical patient data. Each of the 160 patient images, the PET and sMRI scan for each patient, were first manually analysed using version 12 of the Statistical Parameteric Mapping [70] (SPM) software; the origin of the image was set to coincide with the anterior commissure. Next, the sMRI images for each patient were pre-processed for connectome-graph visualization. The SPM software was used, on each patient sMRI, to perform a unified segmentation procedure [70]. The unified segmentation procedure identifies grey matter, white matter, cerebrospinal fluid, skull, and exterior regions. Following this, the grey and white matter segmentations served as input for spatial normalisation using the DARTEL [71] toolbox; the outputs of which were a composite template for the Aβ patient group and a separate composite template for the τP patient group. All patient grey matter images were then normalised to MNI-152 space using their group-specific optimised DARTEL template.

The next step of the pipeline is to treat the PET images for both the Aβ and τP groups. This step relies on the fact that we have already manually relocated the origin of the PET images to the approximate visual location of the anterior commissure as mentioned above. The first new step for this portion of the pipeline is to use the SPM software to co-register the PET images to their sMRI counterparts. This co-registration step is the genesis of the original data procurement requirement that an sMRI scan is conducted no later than one year beyond the PET image acquisition date. The coregistered PET images are then spatially normalised using the DARTEL template and corresponding subject deformation fields derived from the sMRI pipeline (c.f. above). Finally, SUVR values were computed, using SPM, by means of a whole-cerebellar reference region; the skull was then stripped. A voxel-wise mean, across all subjects, was taken to produce a representative SUVR map of both AD cohorts. This completes the first portion of the connectome-graph based imaging analysis pipeline; the result of this step, for both the Aβ and τP group, is shown in the top row of Fig 13.

[Figure omitted. See PDF.]

Fig 13. Skull-stripped, cross-sectional Alzheimer’s patient cohort SUVR intensity.

Top row: averaged SUVR data is shown. Bottom row: top 30% of SUVR intensities are visible. For both rows: (left side) 18F-AV45 florbetapir Aβ radiotracer SUVR and (right side) 18F-AV-1451 flortaucipir τP radiotracer SUVR. Darker colors correspond to higher SUVR values.

https://doi.org/10.1371/journal.pcbi.1008267.g013

The averaged SUVR data of Fig 13 (top row) reports a general view of uptake across the whole brain. In order to visualize significant features of the data: the skull-stripped SUVR image volumes, in NIfTI file format, are visualized using Paraview vis-a-vis the NIfTI Paraview plugin. The volume opacity then set so that the top 30% of the SUVR intensity range in the data is visible; c.f. Fig 13, middle row. Doing so: we immediately see notable features of significance reported in previous radiotracer studies; in particular the familiar [68, 69] temporal and parietal dominance of the τP radiotracer uptake distribution (Fig 13, middle right) are visible. In order to compare simulation results to the patient data of Fig 13 we now employ a connectome-graph data visualization software process. This portion of the general pipeline uses functionality from both SPM and the Nilearn [72] Python library. Regional masks were produced using the Lausanne multiresolution atlas [73] parcellation to the MNI ICBM 152 non-linear 6th generation symmetric volume [74]; generating over 1000 distinct masks. The mask volumes were then applied to isolate the SUVR values for each mask in the parcellation and a regional average SUVR was computed. The computed values were then normalized to lie in the interval [0, 1] by dividing all regional SUVR averages by the global maximum average SUVR intensity; the normalized values, along with the MNI-space coordinates of the region’s centroid were recorded as output. The regional normalised SUVR values and the MNI coordinates were then used as input to the Python application programming interface of the Nilearn [72] software. Using the connectome visualization capabilities of Nilearn we rendered this information using a glass brain view with the highest 30% of values shown; see Fig 14. A comparison with Fig 13 (bottom row) shows that characteristic PET features associated with Alzheimer’s disease [68, 69] are once more prominent in the connectome view of the top 30% of SUVR intensities.

[Figure omitted. See PDF.]

Fig 14. A connectome-graph view of the normalized patient SUVR data.

The (left side) 18F-AV45 florbetapir Aβ radiotracer SUVR and (right side) 18F-AV-1451 flortaucipir τP radiotracer SUVR. Highest 30% of connectome regional values are visible. Darker colors correspond to higher SUVR values.

https://doi.org/10.1371/journal.pcbi.1008267.g014

Demonstrating that the mathematical model of (8)–(11) is capable of achieving distributions of toxic Aβ and τP that resemble the PET data of AD patients is a multi-step process. We note that the demonstration endeavored here is illustrative and does not constitute a full validation of the model; it will, however, fully justify that the fitting of real-world data is within the capacity of the model. First, we set all regions in the connectome to a state of secondary tauopathy with the general synthetic parameters given by those in Table 3. The Aβ–τP interaction parameter, b3, was modified in several regions. All of the modifications to b3 were symmetric; that is, they were made in both the left and right hemispheres of the corresponding region. The modified interaction parameters for connectome vertices in select secondary tauopathy regions are shown in Table 4. Finally, the connectome vertices in a total of five brain regions, in both hemispheres, were put into a state of primary tauopathy by changing the values of b2 and b3 to correspond to states in this regime. The primary tauopathy regions, and their parameters, are listed in Table 5.

[Figure omitted. See PDF.]

Table 3. Comparison with ADNI Alzheimer’s patient PET data.

General Synthetic parameters.

https://doi.org/10.1371/journal.pcbi.1008267.t003

[Figure omitted. See PDF.]

Table 4. Regional interaction parameter variation in secondary tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.t004

[Figure omitted. See PDF.]

Table 5. Primary tauopathy regions and parameters.

https://doi.org/10.1371/journal.pcbi.1008267.t005

The connectome vertex parameters given by Table 3 and regional vertex parameter modifications pursuant to Tables 4 and 5 describe a mixed-modality mathematical model; the connectome graph contains vertices in a state of primary tauopathy and vertices in a state of secondary tauopathy. The model Eqs (8)–(11) were solved with the regional parameters, and modifications, described above. Seeding patterns for both Aβ and τP are identical to those discussed at the beginning of the Alzheimer’s application section. The patient SUVR data, visualised on the conectome, is shown in Fig 14; the results of the simulation are shown in Fig 15 at time t = 78 in accordance with the mean age of the Aβ and τP cross-sectional study parameters (c.f. Table 2); the highest 30% of values are visible. A comparison of Fig 13 (bottom row) and Fig 14 to that of Fig 15 shows that the model can indeed capture salient characteristics of Alzheimer’s disease proteopathy as indicated by SUVR intensity.

[Figure omitted. See PDF.]

Fig 15. Results of a mixed-modality simulation.

(left) Toxic Aβ population and (right) toxic τP population are shown at time t = 78. The top 30% of nodal values are visible; darker colors correspond to higher values.

https://doi.org/10.1371/journal.pcbi.1008267.g015

Thus, this preliminary result clearly demonstrates that the model can recover primary features of Alzheimer’s disease proteopathy and that more mathematically comprehensive analyses are warranted; for instance, investigations using (variational) Bayesian methods [75, 76] may be compelling for further study of the model alongside patient data for cognitively normal, mildly cognitively impaired, early and late onset Alzheimer’s disease cohorts.

A simple model of local and non-local neuronal damage

This section briefly examines the use a simple measure of neuronal damage with the minimal level of complexity necessary to take into account both local and non-local effects. The intent is to explore the qualitative differences between the primary and secondary tauopathy regimes and the effect of varying: the toxification rate of Aβ on τP; and the rate of aggregation due to non-local influence. The continuous Eq (3) were augmented (Methods, A continuous mathematical model) with a coarse-grained damage model (4). We recall that q(x, t) represents a first-order assessment for neuronal cell body damage vis-a-vis a, potentially variegated, set of coupled mechanisms. These mechanisms are not individually differentiated; however, they are assumed to be correlated with the presence of toxic Aβ, toxic τP or with those mechanisms requiring both (c.f. the discussion surrounding Eq (4)).

The damage model (4) has several coefficients: k1 and k2 mediate the damaging effect of toxic Aβ and τP respectively. The rate coefficient k3 reflects damage, such as the rate of neuronal death following over-excitation, resulting from the combined presence of toxic Aβ and toxic τP. Finally, k4 determines the rate of transneuronal damage propagation; thus reflecting aggregate neuronal death as a result of communication disruption to and from regional neighbors.

In this illustrative example we consider the parameters(2)as a baseline from which to begin investigation. These parameters have been chosen to reflect a few clinical observations. First, k1 is chosen as significantly less than k2 to reflect the correlation [9, 10, 12, 13] of toxic τP neurofibrillary tangles with various forms of neuronal damage (e.g. intracellular NFT-induced neuron death, atrophy etc). Second, toxic effects of τP are increased in the presence of toxic Aβ [12, 26, 28, 29, 30, 31, 32, 33] thus, k3 is taken larger than k2.

As a first point of enquiry: we consider our baseline tauopathy patient parameters (Table 1) and vary the deafferentation parameter k4 across three orders of magnitude from the initial value given in (2). Fig 16a and 16b show the results. Note that, in each subfigure, the dashed lines correspond, from left to right, to monotonically decreasing values of k4; the far left dashed curve is k4 = 1.0, the next curve to the right is k4 = 1 × 10−1, the next is k4 = 1 × 10−2, and so forth, down to the final (rightmost) curve corresponding to k4 = 1 × 10−6. In both figures the baseline deafferentation curve, k4 = 1 × 10−3, is instead solid (and red) for emphasis. Fig 16c and 16d show the effect of increasing b3; we have incremented b3 by two, from baseline, for each case. As expected an overall increase in toxic τP, for primary tauopathy and for secondary, is observed with the increase in b3. However, the limiting behavior of the deafferentation baseline coefficient choice, k4 = 1 × 10−3, remains; which justifies our choice of k4 in (2).

[Figure omitted. See PDF.]

Fig 16. Aggregate damage in primary and secondary tauopathy.

Aggregate damage (dashed; except k4 = 1 × 10−3 solid, red) curves in the base primary (a) and secondary (b) tauopathy patients. Damage with increase toxic protein interaction, b3, in primary (c) and secondary (d) tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.g016

The staging of the damage is presented in two figures: primary tauopathy in Fig 17 and secondary tauopathy in Fig 18. Each set of figures includes an overhead horizontal plane view in addition to a sagittal view of the right hemisphere. A visualization starting time was selected to coincide with the first visibility of 5% damage, in any nodes, while an ending time was selected such that the damage progression appeared qualitatively equal. Progression times are uniformly spaced within this interval to allow for a direct comparison between the damage distribution within the two regimes. An immediate observation is that a 5% damage detection is latent within the secondary model, starting at t = 95, compared to the primary tauopathy paradigm at t = 80.

[Figure omitted. See PDF.]

Fig 17. Damage progression in primary tauopathy.

Horizontal plane view (top row) with opacity exaggerated (second row) progression. sagittal view (third row) with opacity exaggerated (fourth row) progression. Dark blue indicates the minimal damage value of q = 0.0; bright red indicates the maximum of q = 1.0. Intermediate values are: purple (q = 0.14), sky blue (q = 0.29), green (q = 0.43), yellow (q = 0.57), orange (q = 0.71), and dark red (q = 0.86).

https://doi.org/10.1371/journal.pcbi.1008267.g017

[Figure omitted. See PDF.]

Fig 18. Damage progression in secondary tauopathy.

Horizontal plane view (top row) with opacity exaggerated (second row) progression. sagittal view (third row) with opacity exaggerated (fourth row) progression. The color scale is identical to that of Fig 17.

https://doi.org/10.1371/journal.pcbi.1008267.g018

It is challenging to discern differences between the fully opaque horizontal views of Fig 17 v.s. Fig 18; some discrepancies are apparent in the sagittal views, however. Relative opacity exaggeration is used to gain further insight. At each time the minimum and maximum damage, denoted Dmin and Dmax, was computed across all regional nodes of the brain connectome; opacity was then set to linearly increase from: fully transparent at the average ; to fully opaque at the maximum value Dmax. The resulting opacity exaggeration scheme shows, at each time step, the relative distribution of the most damaged regions.

The aforementioned opacity scheme leads to a further observations. First, the distribution of relative significant damage in primary tauopathy (Fig 17, second and fourth rows) is clustered more centrally to the toxic τP seeding site of the transentorhinal cortex. Conversely, the distribution of relative significant damage in secondary tauopathy (Fig 18, second and fourth rows) is distributed in the direction of the temporobasal region; a site associated with Aβ seeding. As the disease progresses, t = 103 and t = 114 in Figs 17 and 18 respectively, we see two distinct differences: relative damage is more connected, in the horizontal plane, in addition to more diffuse in the coronal direction, of the sagittal plane, for the case of primary tauopathy; in secondary tauopathy the relative damage in the horizontal plane forms three distinct clusters while severe damage in the sagittal plane is follows the temporobasal and frontomedial directions.

It is increasingly difficult to visually detect qualitative patterns in later stages of significant damage progression; that is, t ≥ 125 for primary tauopathy and t ≥ 133 for secondary. Nevertheless it appears that late stages, t = 148 and t = 170, for primary tauopathy display a more diffuse distribution of significant relative damage away from the transentorhinal region; whereas late secondary tauopathy, t = 151 and t = 170, show more comparative significant damage in the areas associated with Aβ initial seeding.

Taken collectively: these observations suggest that damage onset and the relative distribution of severe damage may offer distinct points of view for application modelling to both typical Alzheimer’s disease along with its neuropathological subtypes [77, 78].

Discussion

In this section we reflect on the analytic and computational results of the manuscript. We first list some advantages and limitations of the current perspective; in we then discuss several results in the context of the current literature and offer further questions and brief concluding remarks. The proposed model is based on physical protein aggregation kinetics; the simplest such two-famiy-two-species interacting protein model one could posit. Nevertheless, the model is mathematically sophisticated enough to evince two distinct pathology regimes, termed primary and secondary tauopathy, of potential clinical interest. After discretizing (3) on a structural connectome: an approachable system of non-linear ordinary differential Eqs (8)–(11), emerges which can be solved using standard mathematical software; such as Mathematica or Matlab. As a result, we expect the model to be widely appealing to the computational neurodegenerative disease community as a starting point for gaining further insight into protein-protein interactions in the context of Alzheimer’s disease.

Advantages and limitations

The deterministic nature of the model, c.f. (3), has at least three distinct advantages over complex, stochastic models: first, the reaction terms of (3) represent simplified [49, 50], but physical, protein aggregation kinetics with a basis in experimental measurement [50, 52, 53, 54, 56]; second, the connectome-discretized Eqs (8)–(11), can also be readily implemented using off-the-shelf mathematical software (e.g. Mathematica or Matlab, etc). Thus, (8)–(11) are easily approachable and do not require probabilistic postulations, based on data or otherwise, regarding underlying distributions. A third advantage is that (8)–(11) are amenable to an a-priori mathematical analysis. This analysis is immutable in nature and much can be observed as a result of using standard methods from the theory of ordinary differential equations and non-linear diffusion-reaction systems. Conversely, probabilistic models may need extensive tuning, reformulation or data curation in order to determine a model’s emergent properties. An independent investigation, i.e. model fitting and application, founded on datasets with differing fidelity may produce divergent results. Such models essentially act in service to deeply mine a set of data but are not always directly helpful to elucidate the impact of individual disease mechanisms.

Conversely, (8)–(11) has inherent limitations. As discussed in recent literature: [35, 79] there are challenges surrounding the acquisition of the parameters in deterministic models such as (8)–(11). In vitro kinetic parameters, regulating the multiplication and growth of several proteins, have been ascertained for Aβ, τP, α-synuclein and others; c.f. the citations in [79]. If we disregard the clearance terms in the prototypical heterodimer model, (8) and (9), then precisely two kinetic coefficients remain: source production (a0) and healthy-to-toxic conversion (a2). The in-vitro experimental estimation of protein-specific aggregation kinetic parameters, however, typically relies on more complex theoretical models: consisting of at least five kinetic parameters; and an infinite number of equations (c.f. for instance [48, Sec. 3.2]). It is therefore not immediately clear how to obtain explicit values for the kinetic rate parameters of (8)–(11). Asymptotic expansions have provided links between the rate coefficients of other more complex models and their simpler counterparts, c.f. [Sec. 2.2] [47], and such approaches may provide insight into reducing experimental parameters from the five-parameter models [79, 48] to those of (8)–(11). A further complication, though, is that even if explicit, experimentally verified, in-vitro parameters were available for (8)–(11) these do not necessarily translate into the correct parameters in vivo [79] where indirect, and locally varying, mechanisms (such as aggregates interacting with cells, or the effects of inflammation on aggregation dynamics) may play a role in altering the associated rates; either globally or regionally.

Observations and open questions

Neurodegenerative diseases are complex and multi-scale processes. The point of view of (8)–(11) is to reduce this complexity by considering a collection of aggregate mechanisms and their implications. For instance (3) can be viewed, more conversationally, as the following collection of general mechanisms: (a) there exists two protein families; (b) each family has a healthy and toxic species; (c) these species are produced and cleared at some aggregate (regional) rate (d) any movement of these species, within the brain, is primarily determined by the macroscale axonal structure; (d) healthy proteins within a family can become toxic, at some (regional) rate, based on the presence of other toxic proteins of that family; and (e) the conversion of healthy-to-toxic proteins, for the second family, is further influenced by the presence of the toxic population of the first family. The current literature suggests that this collection of observations outlines a minimal prion-like model of Alzheimer’s disease progression.

In section the discussion on advantages and limitations it was mentioned that one advantage of simple deterministic models, such as (8)–(11), is that the impact of individual mechanisms can be elucidated and several emergent behaviors can be ascertained a priori. Models such as (3) therefore lead, naturally, to additional questions and serve as a trailhead for further development. The first, and critical, observation is that: (8)–(11) implies that the local balance of clearance, e.g. (23), plays a fundamental role in disease initiation. In light of the seminal work of Braak and Braak [11] this leads naturally to the question: what are the local (toxic τP) clearance properties characterizing the transentorhinal region, (which defines the early Braak stages) and how do these local properties differ from other regions? Aspects of the fine-scale clearance mechanisms of toxic τP remain unclear or are even controversial [80, 81, 82, 83]. Nevertheless, our simple framework reinforces the sentiment echoed by experimentalists: that understanding these processes may be critical to a mechanistic understanding of the initiation of the disease cascade.

A second observation emerging from (8)–(11) is that the progression of Alzheimer’s disease may consist of a confluence of brain regions simultaneously in differing states characterized by contrasting fundamental dynamics. In particular (Methods, Stability and Disease phenomenology) even our simple model of AD development suggests potentially complex disease phenomenology; one where τP can evolve independently of Aβ (termed primary tauopathy) and one where τP depends intrinsically on the presence of Aβ. Furthermore, the line between these two regimes is demarcated by: the balance of local clearance; and the degree of local influence of Aβ on the toxification of τP [12, 26, 30, 31, 32, 33] as expressed by the bulk parameter b3. Depending on these local attributes we could have some areas of the brain in a state of primary tauopathy and others in a state of secondary tauopathy; the latter regions having their tauopathy delayed until a toxic Aβ population is established while the former regions are free to develop toxic τP and NFT independently. This leads naturally to another fundamental line of further enquiry: what are the simplest additional relations, extended (8)–(11), needed to suitably describe the evolution of clearance and toxicity rates alongside protein pathology?

Our simple mathematical model suggests, as a third observation, that the rate of toxic Aβ-τP interaction (i.e. b3) is not a passive facet of disease phenomenology but, rather, may play a much more integral role. We have already discussed, above, that b3 plays a role in secondary tauopathy; it can do this by lowering v4, in (17), thus ensuring that is an admissible state. Other interesting observations regarding b3 were discussed in the results section (c.f. A simplified model of Alzheimer’s disease proteopathy). The observation regarding the impact of b3 in local and non-local neuronal damage formation is straightforward so we will mention it first; namely, that transneuronal damage propagation in the model has a ‘minimum speed’ in both primary and secondary tauopathy; and that increasing b3 has no effect. In particular Fig 16, and the surrounding discussion, suggests that lowering the ‘transmission coefficient’, reflected in k4, below a certain threshold does not lower transneuronal damage; and that increasing b3, for a fixed choice of k4, does not increase the overall propagation of damage. These details lead to a somewhat interesting observation: that the rate of neuronal damage from the structural network topology of the brain may exhibit a baseline, or minimum, value; independent of the details driving local damage (e.g. from local toxicity).

We also observed two features of tauopathy directly related to b3; the time of onset and the τP ‘invasion window’. We recall that the time of onset is defined as the first appearance of toxic τP while the ‘invasion window’ is the timespan starting at a 1% toxic τP concentration and terminating when the asymptotic steady-state value is achieved. In primary tauopathy, disease onset time is virtually unaffected by varying b3 whereas increasing levels of b3 shortens the tauopathy invasion window. In addition, the asymptotic concentration value of toxic τP increases with b3 so that, overall, increasing b3 implies that a more severe tauopathy will develop, faster, at a similar starting point in time (c.f. Fig 8). The picture in secondary tauopathy is different. We see, again, that increasing b3 does increase the severity of the tauopathy (Fig 10b); however, this is where the similarities with primary tauopathy end. First, as b3 increases the time of onset decreases (Fig 10b). Second, the invasion window in secondary tauopathy does not decrease monotonically with decreasing b3 (Fig 10c); rather, we see the invasion window start time and end both decrease, with increasing b3, while the start time decay and end time decay, relative to increasing b3, is different (Fig 11). This is the cause of the initial drop, from b3 = 0.75 to b3 = 3.0, of the invasion window in Fig 10c followed by an increase to a steady invasion window length circa b3 = 12. The observation that increased b3 can decrease the time of onset in secondary tauopathy, which requires the presence of Aβ, while also impacting the invasion window time is reminiscent of the effects associated to the presence of particular Apoliprotein E (APOE) allele configurations. For instance: APOE ϵ4 carriers are more likely to develop AD; toxic Aβ production and deposition is more abundant in APOE ϵ4 carriers; and APOE ϵ4 exacerbates Aβ-related neurotoxicity [84].

In the results subsection ‘a mixed model comparison to Alzheimer’s diseased patient data’ we discussed a mixed-modality instantiation of the model (8)–(11), with some regions in a state of primary tauopathy with all others in a state of secondary tauopathy, using hand-selected synthetic parameters. It was demonstrated that such a mixture of states can reproduce salient features seen in Alzheimer’s disease; in particular, the distinct distribution patterns [68, 69] of 18F-AV-1451 radiotracer are clearly observed. This observation suggests that the model of (8)–(11) is sufficiently rich and implies that the undertaking of a comprehensive data fitting and comparison study is both well warranted and an optimistic endeavor. It is interesting to note, though, that the distinction between primary and secondary tauopathy is not simply one of differently-valued parameters; in particular the two states are differentiated by the balance of clearance inequalities (Methods, Stability). In particular, we have and for primary tauopathy; for secondary tauopathy the latter inequality changes sign to . This observation implies that an arbitrary parameter fitting could produce accurate results, compared to data, while still being questionable since the fitting would imply secondary characteristics regarding regional clearance attributes which may or may not hold. It would therefore be beneficial to carefully consider a data-based measure of regional clearance, for both Aβ and τP, when selecting a data fitting method; possibly incorporated as a constraint or as part of a cost functional. Nevertheless, the dual regimes of primary and secondary tauopathy provide a verdant backdrop for further modelling endeavors; both in terms of fitting to clinical imaging data and to probing and modelling possible ties between b3, APOE configuration, and secondary tauopathy.

In the introduction section we mentioned that the role of Tau in AD formation, and development, is beginning to be recognized as a potentially significant factor. Despite this, open questions about the nature of tau, and tauopathy, in AD remain. For instance: one could argue that healthy τP, being bound to microtubules, should not be diffusing at all. The literature suggests that even healthy tau in healthy neurons exhibit mobility within the cell [85], is secreted into the extracellular space [86, 87], and that extracellular tau is taken up by neighboring neurons [88]; even in the absence of pathology. It is not entirely clear what the correct choice of ρ, in (7), should then be. Despite the literature seeming to suggest that ρ > 0 one could still insist that ρ ≪ 1, or possibly even assert that ρ = 0, should be chosen in (7) for the graph Laplacian of (10). Regarding impacts to the model: this perspective alone would not affect any of the analytic observations (Methods, An Analysis of the continuous model). Indeed, if v is the vector whose jth entry is the healthy tau concentration, vj, in node j then if v is constant, or nearly so, then the graph Laplacian applied to v is zero, or nearly so, regardless of the value of ρ. Thus, since all of the nodes in the computational investigations discussed in the results section had their healthy tau populations set to the same constant value: the effect of any healthy tau diffusion in the simulation results there would be expected to be entirely negligible as well.

The nature of the rates for healthy tau production and clearance, in the literature, are also not fully understood. Indeed, the visual confirmation of the mRNA machinery for localized transcription [89] of tau in axons, and growth cones, is less than two decades old; clarifying important aspects of tau clearance, both healthy and diseased, is an ongoing process [80, 82, 83]. Our results indicate that violation of the balance of clearance inequalities, (23), is fundamental for disease initiation and phenomenology; for instance: if healthy tau were not regenerated, so that b0 = 0, then the regime of ‘primary tauopathy’ (which requires that ) would be an impossibility. This would imply that, in the context of our model, that the development of all tauopathies would require an accompanying amyloidopathy and would seem to preclude those tauopathies which are mostly dominated by toxic τP spreading [8]. It has been observed that: tau is expelled from neurons [86], including healthy ones, on a periodic basis [87]; and that tau plays a role in cell signalling, cell polarity, synaptic plasticity and the regulation of genomic stability [90, 91]. These observations, alongside recent work in adult neurogenesis [92], give good reason to suspect that both b0 > 0 and b1 > 0; at least in the healthy brain and early in disease progression. An open question, though, is how these quantities may change with disease progression. For instance, one could extend the current model by coupling b0 with the damage coefficient q; reflecting the fact decreased healthy tau synthesis could result from neuronal loss, and the decline of neurogenesis [92], throughout AD progression.

The final observation we mention regards an open question surrounding the imaging, and construction, of the structural connectomes used in such network models. We have used an often cited connectome [58, 59]. This connectome is available in various resolutions; the lowest of which consists of 83 vertices (regions of interest) while the highest resolution case, which we have used here, consists of 1015 vertices. However, there are apparent differences in both Aβ and τP staging when solving equations, such as (8)–(11), on the low versus high resolution connectomes. We used the simple, illustrative parameters described for primary tauopathy (Results, A simplified model of Alzheimer’s disease proteopathy) and recorded the average regional tau concentration at six fixed time points; the time points were selected to span disease progression. Fig 19 shows the results of this tau staging experiment for nine regions.

[Figure omitted. See PDF.]

Fig 19. Toxic τP average regional concentration; six fixed time points.

83 (left) versus 1015 (right) vertex connectomes.

https://doi.org/10.1371/journal.pcbi.1008267.g019

In Fig 19, we can see that the two different resolutions of connectomes, derived from the same set of patient data, offer distinct staging patterns for tau progression. This implies that the connectome itself may play a significant role in retrieving results that match clinical data. Since the parameters of (8)–(12) have a physiological interpretation a simple ‘fitting’ to available clinical data is not satisfactory. Fig 19 suggests that developing a more rigorous understanding of computational staging behavior should be endeavored seriously and from first principles. Validating computational (tau) staging behavior, at different connectome resolutions, against clinical standardized-uptake-value-ratio (SUVR) studies, e.g. [93], is an important next step.

Concluding remarks

We have presented a novel, minimal, and deterministic theoretical mathematical model of protein propagation that includes two interacting protein species. The model is motivated by recent experimental evidence regarding the potential importance that interactions between Aβ and τP may play in the development of AD pathology [14, 16, 23, 24, 25]. The primary contributions of the current manuscript are: clearly, and mathematically, establishing the intrinsic dependence of the model on the balance of clearance inequality, (23), and the stability analysis of the modes of primary and secondary tauopathy (Methods, Stability and Disease Phenomenology); and establishing the speed of propagation of toxic fronts (Methods, Front propagation). Further novel contributions of interest include: demonstrating qualitative properties of disease propagation and damage, in primary and secondary tauopathy (Results, A simplified model of Alzheimer’s disease proteopathy), using globally constant, but non-physical, parameters; and demonstrating that the model can achieve tau distributions that reflect canonical patterns in Alzheimer’s patient data (Results, A mixed model comparison to Alzheimer’s diseased patient data). In particular, we have seen that the topology of the brain connectome leads to complex behavior in both pathological regimes. Finally, we have contextualized numerous analytic and computational observations with reference to the current literature and drawn attention to open avenues of further research suggested by the current work.

Alzheimer’s disease is a complex and multi-scale disease. The need for mathematical models, presenting observed disease characteristics, that are computationally tractable is pressing. Our findings suggest that further enquiry into both protein interaction and clearance processes is an important path forward in elucidating key mechanisms in the progression of these diseases. Due to the ease of implementation of (8)–(12), and the widespread interest in computational neurodegenerative disease, we hope that this model will be appealing, to the community, for probing the nuances of protein-protein interactions in neurodegenerative disease development.

Methods

A continuous mathematical model

The simplest possible deterministic aggregation model accounting for the interaction of two protein families, each consisting of a healthy and toxic population, is the heterodimer model [49, 51, 50]. In the heterodimer model: a toxic, misfolded seed protein recruits a healthy protein, induces misfolding, and then fragments; ultimately producing two copies of the misfolded toxic variant. The heterodimer model views these three processes as a single step and expresses this step as an overall mean rate of reaction; such rates can be determined experimentally. Indeed, measuring the mean rates of protein self-aggregation mechanisms, providing best-fit mean aggregation dynamics to deterministic models such as the heterodimer model, is a thriving field of contemporary research [50, 52, 53, 54, 55, 56] and our choice of a deterministic model is inspired by such work. Fig 20 demonstrates the primary molecular mechanism of the heterodimer model; the healthy (blue) protein is approached by the toxic (red) protein and undergoes three separate transitions (small arrows) which are treated as a single transition (long arrow) taking a healthy protein to a misfolded, toxic state.

[Figure omitted. See PDF.]

Fig 20. Kinetics of the heterodimer model.

Healthy protein (blue) and misfolded toxic protein (red) transition to two toxic proteins (long arrow) via, from left to right, the kinetics of: recruitment, induced misfolding, and fragmentation.

https://doi.org/10.1371/journal.pcbi.1008267.g020

We are interested in the interaction between two different protein families; motivated by the Aβ and τ interactions observed in AD. Towards this end we will consider two heterodimer models: one for Aβ and one for τP. These two models will be coupled together by a single term reflecting that the formation of new toxic τP can be enhanced by the presence of Aβ. The heterodimer model was originally posed [49, 50] as a continuous, non-linear, partial differential reaction-diffusion equation for a single protein. To define the model for our two protein families: let be a spatial domain of interest and, for x ∈ Ω and time , we denote by u = u(x, t), and v = v(x, t) the concentration of healthy Aβ and τP. Similarly, we denote by , and , the concentration of toxic Aβ and τP, respectively. Then, the concentration evolution is governed by(3)The first two equations, above, correspond to the usual heterodimer model for the healthy and toxic variants of the protein u; note that the second heterodimer model, for the variants of protein v, deviates from the form of the first by a single balanced term with coefficient b3. The system (3) could apply to any two families of interacting proteins; though the model is inspired by AD. The parameters are as follows: (a0, b0) are the mean production rates of healthy proteins, () are the mean clearance rates of healthy and toxic proteins, and (a2, b2) reflect the mean conversion rates of healthy proteins to toxic proteins. The coupling between the two, otherwise separate, heterodimer models for Aβ and τP, is realized via b3. The b3 predicated terms arise from the mode of interaction assumption, c.f. M1 above, dictating that the presence of Aβ augments the conversion process of healthy τP to toxic τP. We note that toxic Aβ acts as an enzyme in this process and is therefore not depleted. In the absence of production and clearance maps, we assume that all these parameters are constant in space and time. The symmetric diffusion tensors D1,2 and characterize the spreading of each proteins. For isotropic diffusion, these tensors are a multiple of the identity, D1,2 = d1,21 and ∇ ⋅ (D1,2 ⋅ ∇(u)) = d1,2Δ(u) is the usual Laplacian operator (similarly for , v and ) For anisotropic diffusion, the eigenvector with the largest eigenvalue describes the direction of faster diffusion which is used to model preferential propagation along axonal pathways [37].

The coupled system of Eq (3) dictates the spread, genesis, and clearance of two healthy species, u and v, and two toxic species, and , of proteins throughout the domain Ω. The presence of toxic proteins near a point x ∈ Ω can disrupt the extracellular environment of neurons near x and impair their intracellular function. A broad range of coupled effects can contribute to neuronal impairment; including: chronic inflammation, erosion of the blood-brain barrier surrounding vessels, accelerating tau hyperphosphorylation, disrupting normal synaptic efficacy, and deafferentation, among others. A hallmark of neurodegenerative proteopathies is cognitive decline; propelled by the various coupled effects induced by the presence of toxic aggregates and the widespread erosion of neuronal integrity. The nuanced coupling between these disparate deleterious effects is not well understood; nevertheless, we employ the observation that such effects are generally correlated with larger concentrations of misfolded aggregates to define a gross measure of regional neuronal ‘damage’ denote by q(x, t) ∈ [0, 1]. This damage variable takes the perspective that q(x, t) = 0 signifies that the neurons in a neighborhood of x ∈ Ω are functional and healthy whereas q(x, t) = 1 implies that neurons near x have reached a fully-degenerate asymptotic state whereby they are either no longer functioning or fully deceased. For the evolution of the damage we assume a simple, first-order rate model:(4)When k4 = 0: the evolution Eq (4) can be seen as a first order reaction model, i.e. exponential decay, for the transformed variable ; the associated rate of decay is then dependent on the deposition concentration of pattern of and . The first two parameters in (4) denote contributions to neuronal dysfunction, near x, due to presence of isolated toxic aggregates; while the third term accounts for contributions requiring, or accelerated by, the presence of both toxic species aggregate species together. Thus, the third term engenders both toxic effects M2 and M3; c.f. the introduction section. The first three terms of (4) account for neuronal dysfunction in a neighborhood of the point x while the last term, , incorporates non-local contributions, such as transneuronal degeneration, whereby the impairment, or death, of neighboring neurons can increase [57] the probability of impairment near x; thus leading to an increased mean rate of local decline. This nonlocal term does not have a simple representation within the continuous framework as the positions of neuronal bodies is not explicitly encoded. However, we will see that in the discrete case, there is a natural way to take this effect into account and we will delay the discussion of this term until the next section.

A network mathematical model

A simple coarse-grain model of the continuous system can be obtained by building a network from brain data. The construction is obtained by defining nodes of the network to be regions of interest in the domain Ω, typically associated with well-known areas from a brain atlas. The edges of this network represent axonal bundles in white-matter tracts. The brain connectome is then modeled as a weighted graph with V nodes and E edges obtained from diffusion tensor imaging and tractography. A network approximation of the diffusion terms, having the general form ∇ ⋅ (D∇u) or similar, in the system (3) will be constructed by means of a weighted graph Laplacian. The weights of the weighted adjacency matrix W, used to construct the graph Laplacian, are selected as the ratio of mean fiber number nij by mean length squared, , between node i and node j. That is:(5)The choice of weights, above, are consistent with the inverse length-squared dependence incurred by canonical discretizations of the continuous Laplace (diffusion) operator appearing in (3). The weighted degree matrix is the diagonal matrix with elements(6)Additionally, we define the graph Laplacian L as(7)where ρ is an overall effective diffusion constant. The adjacency matrix for the simulation is derived from the tractography of diffusion tensor magnetic resonance images corresponding to 418 healthy subjects of the Human Connectome Project [58] given by Budapest Reference Connectome v3.0 [59]. The graph contains V = 1015 nodes and E = 70,892 edges and is shown in Fig 21.

[Figure omitted. See PDF.]

Fig 21. A high-resolution brain structural connectome graph.

(Bottom left) The average of 419 brain connectomes with V = 1, 015 vertices spanning (bottom right) 49 associated brain regions; the strongest 2,773 edge connections are shown. The weighted adjacency matrix (top) corresponding to the averaged connectome (bottom).

https://doi.org/10.1371/journal.pcbi.1008267.g021

Let be the concentration of healthy and toxic Aβ and denote the concentration of healthy and toxic τP at node j. The network equations corresponding to the continuous model then take the form of a system of first-order ordinary differential equations. There are four such equations, , for each of the 1,015 vertices in the system; these four nodal equations are: (8)(9)(10)(11)where j = 1, …, V = 1, 015. Similarly, for the damage model we define a damage variable qj at each node j and assume the same law(12)where Ajk is the weighted network adjacency matrix Ajk = njk/ljk if j ≠ k (and njk > 0) and 0 otherwise. Thus k4 has the interpretation of a ‘transmission speed’; the time it takes for the effects of degeneracy in cell k to reach cell j. The weighting chosen in the adjacency matrix term is inspired by the propagation of transneuronal degeneration from a node to its neighbors.

An analysis of the continuous model

Homogeneous system.

It is instructive to start with an analysis of the homogeneous system obtained by assuming that there is no spatial dependence. This analysis applies to both network and continuous models. In this case, both systems reduce to the dynamical system(13)where all variables and initial conditions are assumed to be positive and all parameters are strictly positive.

Damage evolution.

For the homogeneous system above the concentrations remain homogeneous for all time. Damage, in contrast, is node-dependent and expressed by the (nodal) variable qj ∈ [0, 1]. Indeed, in this case, the non-local term associated with transneuronal degeneration, commensurate with the tensor Ajk in Eq (12), cannot be homogeneous. Nevertheless, the damage dynamics are simple enough to describe. Damage will initially increases linearly in time, homogeneously, from the initial value qj = 0. The increase will then trend exponentially at each node, with node-dependent time scales depending on the local node’s degree, and saturate to the value qj = 1 asymptotically in time at each node.

Stationary points.

The stationary points and stability of the homogeneous system (13) are instructive; they inform the disease dynamics implied by the local model. The system (13) can exhibit one, two, three, or four stationary points depending on the parameters; these are:

1. Healthy τP-healthy Aβ: This stationary state is always a solution to (13) and is descriptive of an individual with zero toxic load; no amyloid plaques or neurofibrillary tau tangles. The state is given by:(14)

2. Healthy τP–toxic Aβ: This state describes a diseased brain wherein some Aβ plaques exist but the tau fibril (NFT) concentration or that of hyperphosphorylated tau is non-existent or negligible. A description of this stationary state in terms of the base problem parameters is:

In terms of u1 = a0/a1, from (14), and it is given by(15)

Since the concentrations must be non-negative: the form of , above, implies that u1 ≥ u2. This results in the condition of . In other words either the clearance term of toxic Aβ must be sufficiently small, the conversion term must be sufficiently large, or a ratio of the two, to allow for the existence of a toxic state.

3. Toxic τP–healthy Aβ: This stationary state is a conceptual dual to the previous state above; granted, toxic τP does not influence the Aβ population whereas Aβ does induce additional τP formation. As in (15) we express this state, immediately here, in terms of u1 = a0/a1 and v1 = b0/b1 as(16)Requiring v1 ≥ v3 implies that .

4. Toxic τP–toxic Aβ: This stationary state reflects the invasion of a patient’s brain by both toxic amyloid beta and toxic tau. As in (14)–(16) we write the state in terms of the previous state variables u1 = a0/a1, , , and v4, defined below, as:(17)Introducing(18)into (17) gives(19)

Stability.

We briefly discuss the stability of the stationary points. In addition we distinguish between the two possible ‘disease’ phenomena of (13): the case of a disease system characterized by the dynamics of a four-stationary-point model and the case of a disease system characterized by three fixed points.

Eigenvalues of the linearized system.

The linearization of (13) about any fixed point is governed by the Jacobian matrix(20)The first two eigenvalues of (20) correspond to the Aβ subsystem, e.g. , of (13). Since the coupling of (13) is a one-way coupling these eigenvalues are given by the corresponding eigenvalues of the uncoupled heterodimer model:(21)where and . The third and fourth eigenvalues of (20), corresponding to the coupled tau system of (13), can be written as(22)with and . The form of the tau eigenvalues coincides with those for Aβ when b3 = 0 or when vanishes.

Disease phenomenology.

We can interpret the different stationary states in terms of disease dynamics and define, accordingly, different disease states.

The healthy brain.

A healthy patient represents an instantiation of the healthy stationary state whereby . For the Healthy τP-healthy Aβ state to exist we must have a0 ≤ a1 and b0 ≤ b1, i.e., (u, v) ∈ [0, 1] × [0, 1] are valid concentrations. A failure in healthy clearance, either with an amyloid clearance value satisfying 0 ≤ a1 < a0 or with a tau clearance of 0 ≤ b1 < b0, implies the non-existence of a physically relevant healthy state (c.f. (14). It is instructive to note that the expressions a0/a1 and (respectively b0/b1 and ) express a balance of healthy Aβ production to clearance and toxic Aβ production to clearance (respectively healthy τP and toxic τP production to clearance). Consider the following balance of clearance inequalities:(23)A patient satisfying (23) enjoys full stability to perturbations while in the healthy state (14). That is: if (23) holds with given by (14) then the real parts of the eigenvalues (21) and (22) are negative and the production of small amounts of toxic Aβ, or of toxic tau, or the excess production of healthy Aβ, or healthy tau, results in a quick return to the healthy homeostatic baseline state of (14). The above implies that the model (8)–(11) recognizes the critical role that clearance plays in neurodegenerative diseases. A low value of toxic clearance , respectfully , with sustained healthy clearance or a low value of healthy clearance a1, respectfully b1, with sustained toxic clearance is enough to trigger an instability capable of driving the system away from the healthy state.

The susceptible brain.

From the previous discussion, we conclude that an unfavorable alteration in clearance mechanisms not only renders the healthy state unstable to perturbations but brings into existence the other stationary points characterizing various pathological conditions.

Indeed, a well established clinical biomarker for Alzheimer’s disease is a drop in soluble amyloid concentration in the cerebrospinal fluid; directly suggesting a decrease in a1. Recent evidence also suggest [60] that toxic tau filaments in chronic traumatic encephalopathy patients enclose hydrophobic molecules which may contain blood-born pathogens; a possible result of vascular damage from an impact. Such a finding could imply, for instance, that repeated traumatic injury causes vessel rupture and a subsequent proclivity for this unique form of toxic tau production. The stage is then set to trigger a pathological decline when the critical relation (23), corresponding to tau, is violated due to a balance of increased toxic load and age-induced clearance deficit.

The moment of susceptibility occurs when the inequality of (23) becomes an equality. Mathematically, this parameter configuration is a transcritical bifurcation for the homogeneous system (13) at the coincidence of a combination of the states (14)–(16). Clinically, this is the point whereby additional stationary states are physically meaningful and pathology development becomes a possibility.

The proteopathic brain.

The proteopathic brain has suffered a perturbation from the healthy stationary state; due to the instability in the system this patient is progressing towards a diseased state. The potential pathology phenotypes depend on the patient’s individual parameter values. In particular, if holds then the existence of (15) is physically meaningful and if holds then the same is true of (16). It may be the case, depending on the combination of failed clearance subsystems and specific predisposition for toxic loading, that both relations hold simultaneously. A necessary (clinical) existence criterion for the proteopathic stationary point (17) can be observed directly from the equation for in (19): namely(24)This implies that the parameter b3, defining μ in (18), cannot vanish.

Finally since b3 ≠ 0 and the numerator of of v4, in (19), is always non-negative we see that (17) always exists when u1 > u2 and when both v3, v1 ≥ v4 or when both v3, v1 ≤ v4. An important observation is that, though the modeling of the pathology of (17) is tied to that of (15) it is not inextricably tied to (16); this is due to the fact that we may always choose b3, c.f. (18), such v4 is smaller than both v3 and v1. Thus, with a suitably strong Aβ tau-toxification interaction the state (16) is not needed in order to produce tau proteopathy; that is, the model admits a pathology whereby toxic tau is created solely by the presence of toxic Aβ. Therefore, there are two clinically interesting patient proteopathies for our analysis: the case where the patient model consists of all four disease state equilibria, (14)–(17), and the case where the patient model has the three equilibria (14), (15) and (17).

Primary tauopathy.

In this case, all four equilibria exist which requires both and . An example of this dynamic is shown in Fig 22. We see that the presence of toxic Aβ always implies a higher level of τ P. Indeed, we have(25)We refer to this case as primary tauopathy as the invasion due to τP exists independently of Aβ. The effect of Aβ is to increase the concentration of toxic τP and, possibly, increase the associated damage.

[Figure omitted. See PDF.]

Fig 22. Patient pathology dynamics in primary tauopathy.

(Left) Phase plane with four equilibria. Homogeneous dynamics of the toxic states. Note that this is a two-dimensional slice of the four-dimensional phase space. (Right) When four different states co-exist, only the fully toxic state is stable as shown by the time-dynamics plot. (Parameters: a0 = b0 = a1 = a2 = b1 = b2 = 1, , b3 = 1/2).

https://doi.org/10.1371/journal.pcbi.1008267.g022

Secondary tauopathy.

In secondary tauopathy the evolution of τP depends on the primary invasion of Aβ. Parameters corresponding to secondary tauopathy can be obtained by choosing and (hence, ) while taking b3 large enough so that .

It is useful to outline key observations, explored further in the manuscript, regarding the nature of the dependence of τP pathology () on the presence of Aβ pathology () in this regime. We will see (Front propagation, secondary tauopathy) that the onset of toxic τP follows from the presence of toxic Aβ. We will also see that the speed of toxic τP propagation appears to be limited by the speed of propagation of toxic Aβ. This is distinctly different than the case of primary tauopathy where toxic τP and toxic Aβ can evolve separately. Indeed, Fig 22 (left) shows a stationary point of the form and the fully invaded asymptotic states, Fig 22 (right), satisfy ; clearly indicating that the additional coupling of to , in (11) does not limit tau pathology expression, in primary tauopathy, to that of Aβ pathology. In light of the apparent dependence of toxic τP spreading on toxic Aβ propagation in secondary tauopathy it is instructive to enquire whether the asymptotic level of toxic τP pathology concentration be limited by the asymptotic concentration of toxic Aβ? Fig 23 (right) seems to indicate that this is the case.

[Figure omitted. See PDF.]

Fig 23. Patient pathology dynamics in secondary tauopathy.

(Left) Phase plane with three equilibria. (Right) When three different states co-exist, only the fully toxic state is stable as shown by the time-dynamics plot. (Parameters: a0 = b0 = a1 = a2 = b1 = b2 = 1, , , b3 = 3). Note that trajectories are initialized by taking the initial condition ϵ = 0.005 away from an equilibrium point.

https://doi.org/10.1371/journal.pcbi.1008267.g023

However, as discussed in the results section (c.f. secondary tauopathy in ‘A simplified model of Alzheimer’s disease proteopathy’), the asymptotic state value of is not limited by that of . Fig 10a shows that is possible depending on the parameter selection. This observation leads naturally to a further question: is it a strict requirement of secondary tauopathy to have a toxic Aβ concentration at all? That is: can we have ? Suppose this is the case: then, according to (14) and (17), we have so that u1 = u2. When u1 = u2 we have, again from (17), that v4 = v3; this leaves , in (17), in an indeterminate form. Thus, to understand when u4 = 0, we use that , c.f. (14), and consider:Due to the indeterminate form of the limit above: we apply L’Hôspital’s rule and compute, after simplifying, that(26)As discussed above: the regime of secondary tauopathy occurs when ; resulting in (16) being an invalid steady state. This directly implies that v1 − v3 < 0 so that the expression (26) is negative. The above shows that if the asymptotic toxic Aβ concentration, , vanishes in the case of secondary tauopathy then the asymptotic toxic tau concentration, is necessarily negative; this cannot occur in a physical system. To summarize: the sustained presence of a toxic τP population, in secondary tauopathy, requires the presence of toxic Aβ; this implies that the computational observations, for secondary tauopathy, regarding the reliance of toxic τP development and perpetuation, on the presence of Aβ are strict and not merely one possible method of toxic τP development.

Front propagation.

We can explore the spatio-temporal behavior of the system by first considering a reduction to one dimension () and subsequently analyzing the spread of toxic protein via the study of traveling waves. From the theory of nonlinear parabolic partial differential equations, we expect pulled fronts that connect one equilibrium state to a different homogeneous state [49].

First, consider the two uncoupled fronts emanating from the healthy state and connecting either to or . To obtain these fronts, we linearize (3) around the healthy state and obtain the decoupled system(27)(28)Starting with initial positive data, the system will develop fronts and the asymptotic selected speed is the minimum possible speed for this linear system [61, 62]. Traveling wave solutions to (27) and (28) are obtained explicitly by first performing a traveling wave reduction (u(x, t) → u(z) with z = x − ct and so on for the other variables) and then looking for linear solutions of the form u = C exp(λz) which leads to a family of possible solution with speeds c = c(λ). The smallest such speed is the selected speed for the asymptotic dynamics. In our case, the front speeds are(29)where and denote the speeds of the front from state i to state j (whenever such a front exists) for the Aβ fields and τP fields , respectively. The front speeds for the second transition are(30)Similarly, if both fields are seeded initially, we have(31)We see that these fronts only exist if and/or which are the conditions for the existence of toxic states found in the previous section. Trivially, a front between two states can only develop if such states exist.

Second, we consider the possibility of fronts propagating from equilibrium state 2 to state 4. To do so, we linearize the equations around and repeat the previous steps to find(32)

Primary tauopathy.

As an example of the interactions between the two fronts, we consider a toxic Aβ front on the real axis x propagating to the right interacting with a τP front propagating to the left (see Fig 24). They evolve initially with constant speeds and respectively (Fig 24 top). However, when they overlap, the interaction creates an increase in the concentration of τP (Fig 24 top) which both boosts the front to speed and initiates a new front propagating backward to fill the interval to the global stable equilibrium with speed . The Aβ front is never affected by the presence of toxic τP.

[Figure omitted. See PDF.]

Fig 24. Front dynamics, primary tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.g024

In Fig 24, an Aβ front propagating to the right travels towards a τP front propagating to the left (Fig 24a). The interaction (Fig 24b) increases the toxic level of τP and creates a second front propagating to the right connecting to (Fig 24c). The front profiles are shown at time t = 30, 200, 250, 310. Parameters are as in Fig 22: a0 = b0 = a1 = a2 = b1 = b2 = 1, , b3 = 1/2, which leads to and . Neumann boundary conditions were used on both sides of the finite interval for all variables.

Secondary tauopathy.

As a second example, we consider the case where the Aβ front causes the creation of a non-zero toxic τP state (see Fig 25). Fig 25 depicts the toxic front dynamics of and ; an Aβ front () propagating to the right in a domain with negligible toxic τP (but with a healthy τ population). The passage of the front leads to the rapid expansion of toxic τP (Fig 25b) which evolves at a speed close to (Fig 25c). Hence, it eventually catches up with the front (d) and matches its speed. The front profiles are shown at time t = 180, 200, 215, 300. Parameters as in Fig 23: a0 = b0 = a1 = a2 = b1 = b2 = 1, , , b3 = 3, which leads to and . Homogeneous Neumann boundary conditions are used on both sides of the finite interval for all variables. Initial seeding of toxic τP on the positive interval only with for x > 0 and 0 otherwise.

[Figure omitted. See PDF.]

Fig 25. Front dynamics, secondary tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.g025

Finally, the front propagating from equilibrium state 3 to state 4 is constrained by the evolution of the u and fields. Therefore, we find(33)

Supporting information

S1 Video. Front dynamics in primary tauopathy, synthetic.

Visualization of an illustrative primary tauopathy model simulation using a synthetic channel domain.

https://doi.org/10.1371/journal.pcbi.1008267.s001

(MP4)

S2 Video. Front dynamics in primary tauopathy, connectome.

Visualization of an illustrative primary tauopathy model simulation using a high-resolution structural brain connectome domain.

https://doi.org/10.1371/journal.pcbi.1008267.s002

(MP4)

S3 Video. Front dynamics in secondary tauopathy, synthetic.

Visualization of an illustrative secondary tauopathy model simulation using a synthetic channel domain.

https://doi.org/10.1371/journal.pcbi.1008267.s003

(MP4)

S4 Video. Front dynamics in secondary tauopathy, connectome.

Visualization of an illustrative secondary tauopathy model simulation using a high-resolution structural brain connectome domain.

https://doi.org/10.1371/journal.pcbi.1008267.s004

(MP4)

S5 Video. Toxic Aβ in primary and secondary tauopathy.

Visualization of the toxic Aβ species in an illustrative model simulation using a high-resolution structural brain connectome domain. Qualitative propagation of the Aβ species was identical in both primary and secondary tauopathy.

https://doi.org/10.1371/journal.pcbi.1008267.s005

(MP4)

S6 Video. Toxic τP in primary tauopathy, connectome.

Visualization of the toxic τP species in an illustrative model simulation of primary tauopathy using a high-resolution structural brain connectome domain.

https://doi.org/10.1371/journal.pcbi.1008267.s006

(MP4)

S7 Video. Toxic τP in secondary tauopathy, connectome.

Visualization of the toxic τP species in an illustrative model simulation of secondary tauopathy using a high-resolution structural brain connectome domain.

https://doi.org/10.1371/journal.pcbi.1008267.s007

(MP4)

S1 Data. Structural brain connectome, low-resolution.

A low-resolution graph of the structural brain connectome. The graph is expressed in a standard format (graphml) based on the human-readable XML markup language. This graph consists of 83 vertices (anatomical regions of interest) and 1,654 edges; data for the graph was sourced from freely-available patient connectome data (https://braingraph.org).

https://doi.org/10.1371/journal.pcbi.1008267.s008

(GRAPHML)

S2 Data. Structural brain connectome, high-resolution.

A high-resolution graph of the structural brain connectome. The graph is expressed in a standard format (graphml) based on the human-readable XML markup language. This graph consists of 1,015 vertices (anatomical regions of interest) and 70,892 edges; data for the graph was sourced from freely-available patient connectome data (https://braingraph.org).

https://doi.org/10.1371/journal.pcbi.1008267.s009

(GRAPHML)

S1 Appendix. Additional detail and solver verification.

This appendix contains additional detail regarding the expected behaviour of the network dynamics for specific choices of synthetic parameters. This appendix refers to the equations, and figures, in the main manuscript and is intended to provide a set of verification cases for use in implementing solvers.

https://doi.org/10.1371/journal.pcbi.1008267.s010

(PDF)

Acknowledgments

Data used in the preparation of this article were obtained from the Alzheimer’s Disease Neuroimaging Initiative (ADNI) database (adni.loni.usc.edu). The ADNI was launched in 2003 as a public-private partnership, led by Principal Investigator Michael W. Weiner, MD. The primary goal of ADNI has been to test whether serial magnetic resonance imaging (MRI), positron emission tomography (PET), other biological markers, and clinical and neuropsychological assessment can be combined to measure the progression of mild cognitive impairment (MCI) and early Alzheimer’s disease (AD). For up-to-date information, see www.adni-info.org. A complete listing of ADNI investigators can be found at: http://adni.loni.usc.edu/wp-content/uploads/how_to_apply/ADNI_Acknowledgement_List.pdf.

Citation: Thompson TB, Chaggar P, Kuhl E, Goriely A, for the Alzheimer’s Disease Neuroimaging Initiative (2020) Protein-protein interactions in neurodegenerative diseases: A conspiracy theory. PLoS Comput Biol 16(10): e1008267. https://doi.org/10.1371/journal.pcbi.1008267

References

1. Alzheimer A. Uber eine eigenartige Erkrankung der Hirnrinde. Zentralbl Nervenh Psych. 1907;18:177–179.

2. Stelzmann RA, Norman Schnitzlein H, Reed Murtagh F. An English translation of Alzheimer’s 1907 paper,“Über eine eigenartige Erkankung der Hirnrinde”. Clinical Anatomy: The Official Journal of the American Association of Clinical Anatomists and the British Association of Clinical Anatomists. 1995;8(6):429–431.

3. Walker LC, Jucker M. Neurodegenerative diseases: expanding the prion concept. Annual review of neuroscience. 2015;38:87–103. pmid:25840008

4. Goedert M, Masuda-Suzukake M, Falcon B. Like prions: the propagation of aggregated tau and α-synuclein in neurodegeneration. Brain. 2017;140(2):266–278. pmid:27658420

5. Hardy JA, Higgins GA. Alzheimer’s disease: the amyloid cascade hypothesis. Science. 1992;256(5054):184–186. pmid:1566067

6. Hardy J, Allsop D. Amyloid deposition as the central event in the aetiology of Alzheimer’s disease. Trends in pharmacological sciences. 1991;12:383–388.

7. Selkoe DJ, Hardy J. The amyloid hypothesis of Alzheimer’s disease at 25 years. EMBO molecular medicine. 2016;8(6):595–608. pmid:27025652

8. Götz J, Halliday G, Nisbet RM. Molecular pathogenesis of the tauopathies. Annual Review of Pathology: Mechanisms of Disease. 2019;14:239–261.

9. Cho H, Choi JY, Hwang MS, Kim YJ, Lee HM, Lee HS, et al. In vivo cortical spreading pattern of tau and amyloid in the Alzheimer disease spectrum. Annals of neurology. 2016;80(2):247–258. pmid:27323247

10. Jack CR, Bennett DA, Blennow K, Carrillo MC, Dunn B, Haeberlein SB, et al. NIA-AA Research Framework: Toward a biological definition of Alzheimer’s disease. Alzheimer’s & Dementia. 2018;14(4):535–562.

11. Braak H, Braak E. Neuropathological stageing of Alzheimer-related changes. Acta neuropathologica. 1991;82(4):239–259. pmid:1759558

12. DeVos SL, Corjuc BT, Commins C, Dujardin S, Bannon RN, Corjuc D, et al. Tau reduction in the presence of amyloid-β prevents tau pathology and neuronal death in vivo. Brain. 2018;141(7):2194–2212. pmid:29733334

13. Busche MA, Wegmann S, Dujardin S, Commins C, Schiantarelli J, Klickstein N, et al. Tau impairs neural circuits, dominating amyloid-β effects, in Alzheimer models in vivo. Threshold. 2019;30(40):50.

14. Walker LC, Lynn DG, Chernoff YO. A standard model of Alzheimer’s disease? Prion. 2018;12(5-6):261–265. pmid:30220236

15. Prusiner SB. Prions. Proceedings of the National Academy of Sciences. 1998;95(23):13363–13383.

16. Jucker M, Walker LC. Propagation and spread of pathogenic protein assemblies in neurodegenerative diseases. Nature neuroscience. 2018;21(10):1341. pmid:30258241

17. Olsson TT, Klementieva O, Gouras GK. Prion-like seeding and nucleation of intracellular amyloid-β. Neurobiology of disease. 2018;113:1–10. pmid:29414379

18. Jucker M, Walker LC. Self-propagation of pathogenic protein aggregates in neurodegenerative diseases. Nature. 2013;501(7465):45–51. pmid:24005412

19. Clavaguera F, Lavenir I, Falcon B, Frank S, Goedert M, Tolnay M. “Prion-like” templated misfolding in tauopathies. Brain Pathology. 2013;23(3):342–349. pmid:23587140

20. Goedert M. Alzheimer’s and Parkinson’s diseases: The prion concept in relation to assembled Aβ, tau, and α-synuclein. Science. 2015;349(6248):1255555.

21. Mudher A, Colin M, Dujardin S, Medina M, Dewachter I, Naini SMA, et al. What is the evidence that tau pathology spreads through prion-like propagation? Acta neuropathologica communications. 2017;5(1):99. pmid:29258615

22. DeVos SL, Corjuc BT, Oakley DH, Nobuhara CK, Bannon RN, Chase A, et al. Synaptic tau seeding precedes tau pathology in human Alzheimer’s disease brain. Frontiers in neuroscience. 2018;12:267. pmid:29740275

23. Ittner LM, Götz J. Amyloid-β and tau—a toxic pas de deux in Alzheimer’s disease. Nature Reviews Neuroscience. 2011;12(2):67.

24. Jack CR Jr, Knopman DS, Jagust WJ, Petersen RC, Weiner MW, Aisen PS, et al. Tracking pathophysiological processes in Alzheimer’s disease: an updated hypothetical model of dynamic biomarkers. The Lancet Neurology. 2013;12(2):207–216.

25. Kara E, Marks JD, Aguzzi A. Toxic protein spread in neurodegeneration: reality versus fantasy. Trends in molecular medicine. 2018.

26. Bennett RE, DeVos SL, Dujardin S, Corjuc B, Gor R, Gonzalez J, et al. Enhanced tau aggregation in the presence of amyloid β. The American journal of pathology. 2017;187(7):1601–1612. pmid:28500862

27. Pooler AM, Hyman BTea. Amyloid accelerates tau propagation and toxicity in a model of early Alzheimer’s disease. Acta Neuropathologica Communications. 2015;3(1):14. pmid:25853174

28. Bennett RE, Robbins AB, Hu M, Cao X, Betensky RA, Clark T, et al. Tau induces blood vessel abnormalities and angiogenesis-related gene expression in P301L transgenic mice and human Alzheimer’s disease. Proceedings of the National Academy of Sciences. 2018;115(6):E1289–E1298.

29. Laurent C, Buée L, Blum D. Tau and neuroinflammation: what impact for Alzheimer’s disease and tauopathies? Biomedical journal. 2018;41(1):21–33. pmid:29673549

30. Small SA, Duff K. Linking Aβ and tau in late-onset Alzheimer’s disease: a dual pathway hypothesis. Neuron. 2008;60(4):534–542. pmid:19038212

31. Lloret A, Badia MC, Giraldo E, Ermak G, Alonso MD, Pallardó FV, et al. Amyloid-β toxicity and tau hyperphosphorylation are linked via RCAN1 in Alzheimer’s disease. Journal of Alzheimer’s Disease. 2011;27(4):701–709. pmid:21876249

32. Giraldo E, Lloret A, Fuchsberger T, Viña J. Aβ and tau toxicities in Alzheimer’s are linked via oxidative stress-induced p38 activation: protective role of vitamin E. Redox biology. 2014;2:873–877. pmid:25061569

33. He Z, Guo JL, McBride JD, Narasimhan S, Kim H, Changolkar L, et al. Amyloid-β plaques enhance Alzheimer’s brain tau-seeded pathologies by facilitating neuritic plaque tau aggregation. Nature medicine. 2018;24(1):29. pmid:29200205

34. Walker LC, Jucker M. The exceptional vulnerability of humans to Alzheimer’s disease. Trends in molecular medicine. 2017;23(6):534–545. pmid:28483344

35. Carbonell F, Iturria-Medina Y, Evans AC. Mathematical modeling of protein misfolding mechanisms in neurological diseases: a historical overview. Frontiers in Neurology. 2018;9:37. pmid:29456521

36. Bertsch M, Franchi B, Marcello N, Tesi MC, Tosin A. Alzheimer’s disease: a mathematical model for onset and progression. Mathematical Medicine and Biology. 2016; p. dqw003.

37. Weickenmeier J, Kuhl E, Goriely A. The multiphysics of prion-like diseases: progression and atrophy. Phys Rev Lett. 2018;121(158101).

38. Weickenmeier J, Jucker M, Goriely A, Kuhl E. A physics-based model explains the prion-like features of neurodegeneration in Alzheimer’s disease, Parkinson’s disease, and amyotrophic lateral sclerosis. Journal of the Mechanics and Physics of Solids. 2019;124:264–281.

39. Raj A, Kuceyeski A, Weiner M. A network diffusion model of disease progression in dementia. Neuron. 2012;73(6):1204–1215. pmid:22445347

40. Raj A, LoCastro E, Kuceyeski A, Tosun D, Relkin N, Weiner M, et al. Network diffusion model of progression predicts longitudinal patterns of atrophy and metabolism in Alzheimer’s disease. Cell reports. 2015;10(3):359–369. pmid:25600871

41. Abdelnour F, Voss HU, Raj A. Network diffusion accurately models the relationship between structural and functional brain connectivity networks. Neuroimage. 2014;90:335–347. pmid:24384152

42. Pandya S, Mezias C, Raj A. Predictive model of spread of progressive supranuclear palsy using directional network diffusion. Frontiers in neurology. 2017;8:692. pmid:29312121

43. Poudel GR, Harding IH, Egan GF, Georgiou-Karistianis N. Network spread determines severity of degeneration and disconnection in Huntington’s disease. Human brain mapping. 2019.

44. Iturria-Medina Y, Sotero RC, Toussaint PJ, Evans AC, Initiative ADN. Epidemic spreading model to characterize misfolded proteins propagation in aging and associated neurodegenerative disorders. PLoS computational biology. 2014;10(11):e1003956. pmid:25412207

45. Henderson MX, Cornblath EJ, Darwich A, Zhang B, Brown H, Gathagan RJ, et al. Spread of α-synuclein pathology through the brain connectome is modulated by selective vulnerability and predicted by network analysis. Nature neuroscience. 2019;22(8):1248. pmid:31346295

46. Zheng YQ, Zhang Y, Yau YH, Zeighami Y, Larcher K, Misic B, et al. Connectome architecture, gene expression and functional co-activation shape the propagation of misfolded proteins in neurodegenerative disease. bioRxiv. 2018; p. 449199.

47. S F, A S, A G, Kuhl E. Prion-like spreading of Alzheimer’s disease within the brain’s connectome. Interface R Society. 2019.

48. S F, A S, A G, Kuhl E. Spatially-extended nucleation-aggregation-fragmentation models for the dynamics of prion-like neurodegenerative protein-spreading in the brain and its connectome. J Theor Biol. 2019.

49. Bressloff PC. Waves in neural media. Lecture Notes on Mathematical Modelling in the Life Sciences. 2014.

50. Matthäus F. A comparison of modeling approaches for the spread of prion diseases in the brain. In: Modelling Dynamics in Processes and Systems. Springer; 2009. p. 109–117.

51. F M. The spread of prion diseases in the brain–models of reaction and transport on networks. J Biol Syst. 2009;17:623–641.

52. Cohen SI, Linse S, Luheshi LM, Hellstrand E, White DA, Rajah L, et al. Proliferation of amyloid-β42 aggregates occurs through a secondary nucleation mechanism. Proceedings of the National Academy of Sciences. 2013;110(24):9758–9763.

53. Kundel F, Hong L, Falcon B, McEwan WA, Michaels TC, Meisl G, et al. Measurement of tau filament fragmentation provides insights into prion-like spreading. ACS chemical neuroscience. 2018;9(6):1276–1282. pmid:29590529

54. Meisl G, Yang X, Hellstrand E, Frohm B, Kirkegaard J, Cohen S, et al. Differences in nucleation behavior underlie the contrasting aggregation kinetics of the Aβ40 and Aβ42 peptides. Proceedings of the National Academy of Sciences. 2014;111(26):9384–9389.

55. Meisl G, Kirkegaard J, Arosio P, Michaels T, Vendruscolo M, Dobson C, et al. Molecular mechanisms of protein aggregation from global fitting of kinetic models. Nat Protoc. 2016;11:252–272. pmid:26741409

56. Frankel R, Törnquist M, Meisl G, Hansson O, Andreasson U, Zetterberg H, et al. Autocatalytic amplification of Alzheimer-associated Aβ42 peptide aggregation in human cerebrospinal fluid. Communications biology. 2019;2(1):1–11.

57. Zheng YQ, Zhang Y, Yau Y, Zeighami Y, Larcher K, Misic B, et al. Local vulnerability and global connectivity jointly shape neurodegenerative disease propagation. bioRxiv. 2019; p. 449199.

58. McNab JA, Edlow BL, Witzel T, Huang SY, Bhat H, Heberlein K, et al. The Human Connectome Project and beyond: initial applications of 300 mT/m gradients. Neuroimage. 2013;80:234–245. pmid:23711537

59. Szalkai B, Kerepesi C, Varga B, Grolmusz V. Parameterizable consensus connectomes from the human connectome project: The budapest reference connectome server v3. 0. Cognitive neurodynamics. 2017;11(1):113–116. pmid:28174617

60. Falcon B, Zivanov J, Zhang Wea. Novel tau filament fold in chronic traumatic encephalopathy encloses hydrophobic molecules. Nature. 2019;568:420–423. pmid:30894745

61. Goriely A. A simple solution to the nonlinear front problem. Phys Rev Lett. 1995;75:2047–2050.

62. Van Saarloos W. Front propagation into unstable states. Physics reports. 2003;386(2-6):29–222.

63. Hindmarsh AC, Brown PN, Grant KE, Lee SL, Serban R, Shumaker DE, et al. SUNDIALS: Suite of nonlinear and differential/algebraic equation solvers. ACM Transactions on Mathematical Software (TOMS). 2005;31(3):363–396.

64. Davis TA, Natarajan EP. Algorithm 907: KLU, A Direct Sparse Solver for Circuit Simulation Problems. ACM Trans Math Softw. 2010;37:36:1–36:17.

65. Davis TA. Algorithm 915, SuiteSparseQR: Multifrontal multithreaded rank-revealing sparse QR factorization. ACM Trans Math Softw. 2011;38:8:1–8:22.

66. Grothe MJ, Barthel H, Sepulcre J, Dyrba M, Sabri O, Teipel SJ, et al. In vivo staging of regional amyloid deposition. Neurology. 2017;89(20):2031–2038. pmid:29046362

67. Jucker M, Walker LC. Pathogenic protein seeding in Alzheimer disease and other neurodegenerative disorders. Annals of Neurology. 2011;70:532–540. pmid:22028219

68. Ossenkoppele R, Rabinovici GD, Smith R, Miller BLea. Discriminative Accuracy of [18F]flortaucipir Positron Emission Tomography for Alzheimer Disease vs Other Neurodegenerative Disorders. JAMA. 2018;320(11):1151–1162. pmid:30326496

69. Okamura N, Harada R, Furumoto Sea. Tau PET Imaging in Alzheimer’s Disease. Curr Neurol Neurosci Rep. 2014;14(500). pmid:25239654

70. Ashburner J, Friston KJ. Unified segmentation. NeuroImage. 2005;26(3):839–851. pmid:15955494

71. Ashburner J. A fast diffeomorphic image registration algorithm. NeuroImage. 2007;38(1):95–113. pmid:17761438

72. Abraham A, Pedregosa F, Eickenberg M, Gervais P, Mueller A, Kossaifi J, et al. Machine learning for neuroimaging with scikit-learn. Frontiers in Neuroinformatics. 2014;8:14. pmid:24600388

73. Daducci A, Gerhard S, Griffa A, Lemkaddem A, Cammoun L, Gigandet X, et al. The Connectome Mapper: An Open-Source Processing Pipeline to Map Connectomes with MRI. PLOS ONE. 2012;7(12):1–9.

74. Grabner G, Janke A, Budge M, Smith D, Pruessner J, Collins D. Symmetric Atlasing and Model Based Segmentation: An Application to the Hippocampus in Older Adults. In: Larsen R, Nielsen M, Sporring J, editors. Medical Image Computing and Computer-Assisted Intervention—MICCAI 2006. Berlin, Heidelberg: Springer Berlin Heidelberg; 2006. p. 58–66.

75. Friston KJ, Harrison L, Penny W. Dynamic causal modelling. Neuroimage. 2003;19(4):1273–1302. pmid:12948688

76. Friston KJ, Preller KH, Mathys C, Cagnan H, Heinzle J, Razi A, et al. Dynamic causal modelling revisited. NeuroImage. 2017;199. pmid:28219774

77. Ferreira D, Westman ea Eric. Distinct subtypes of Alzheimer’s disease based on patterns of brain atrophy: longitudinal trajectories and clinical applications. Scientific Reports. 2017;7(1):2045–2322.

78. Jellinger KA. Neuropathological subtypes of Alzheimer’s disease. Acta Neuropathologica. 2012;123(1):153–154. pmid:22009303

79. Meisl G, Knowles T, Klenerman D. The molecular processes underpinning prion-like spreading and seed amplification in protein aggregation. Current Opinion in Neurobiology. 2020;61:58–64. pmid:32092527

80. Tarasoff-Conway J, Carare R, de Leon MJea. Clearance systems in the brain–implications for Alzheimer disease. Nat Rev Neurol. 2015;11(8):457–470. pmid:26195256

81. Iliff JJ, Wang M, Liao Y, Plogg BA, Peng W, Gundersen GA, et al. A paravascular pathway facilitates CSF flow through the brain parenchyma and the clearance of interstitial solutes, including amyloid β. Science translational medicine. 2012;4(147):147ra111–147ra111. pmid:22896675

82. Abbott NJ, Pizzo ME, Preston JE, Janigro D, Thorne RG. The role of brain barriers in fluid movement in the CNS: is there a ‘glymphatic’ system? Acta Neuropathologica. 2019;136(4):388–407.

83. Benveniste H, Liu X, Koundal S, Sanggaard S, Lee H, Wardlaw J. The Glymphatic System and Waste Clearance with Brain Aging: A Review. Gerontology. 2019;65:106–119. pmid:29996134

84. C LC, Kanekiyo T, Xu H, Bu G. Apolipoprotein E and Alzheimer disease: risk, mechanisms, and therapy. Nat Rev Neurol. 2013;2(9):106–118.

85. Konzak S, Thies E, Marx A, Mandelkow EM, Mandelkow E. Swimming against the Tide: Mobility of the Microtubule-Associated Protein Tau in Neurons. J Neurosci. 2017;27(37):9916–9927.

86. Pooler A, Phillips E, Lau D, Noble W, Hanger D. Physiological release of endogenous tau is stimulated by neuronal activity. Embo Reports. 2013;14(4). pmid:23412472

87. Sato C, Barthelemy N, Bateman Rea. Tau Kinetics in Neurons and the Human Central Nervous System. Neuron. 2018;97(6):P1284–1298.

88. Pooler A, Polydoro M, Wegmann S, Nicholls S, Spires-Jones T, Hyman B. Propagation of tau pathology in Alzheimer’s disease: identification of novel therapeutic targets. Alzheimers Res Ther. 2013;5(49). pmid:24152385

89. Aronov S, Aranda G, Behar L, Ginzburg I. Visualization of translated tau protein in the axons of neuronal P19 cells and characterization of tau RNP granules. J Cell Sci. 2002;115:3817–3827. pmid:12235292

90. Morita J, Sobue K. Specification of Neuronal Polarity Regulated by Local Translation of CRMP2 and Tau via the mTOR-p70S6K Pathway. J Biol Chem. 2009;284(40):27734–27745. pmid:19648118

91. Guo T, Noble W, Hanger D. Roles of Tau protein in health and disease. Acta Neuropathol. 2017;133:665–704. pmid:28386764

92. Moreno-Jimenez E, Flor-Garcia M, Llorens-Martin Mea. Adult Hippocampal Neurogenesis Is Abundant in Neurologically Healthy Subjects and Drops Sharply in Patients With Alzheimer’s Disease. Nat Med. 2019;25(4):554–560. pmid:30911133

93. Insel P, Mormino E, Aisen P, Thompson W, Donahue M. Neuroanatomical spread of amyloid β and tau in Alzheimer’s disease: implications for primary prevention. Brain Commun. 2020;2(1):1–11.