Linking Molecular Pathways and Large-Scale Computational Modeling to Assess Candidate Disease Mechanisms and Pharmacodynamics in Alzheimer's Disease

Introduction: While the prevalence of neurodegenerative diseases associated with dementia such as Alzheimer's disease (AD) increases, our knowledge on the underlying mechanisms, outcome predictors, or therapeutic targets is limited. In this work, we demonstrate how computational multi-scale brain modeling links phenomena of different scales and therefore identifies potential disease mechanisms leading the way to improved diagnostics and treatment.

Methods: The Virtual Brain (TVB; thevirtualbrain.org) neuroinformatics platform allows standardized large-scale structural connectivity-based simulations of whole brain dynamics. We provide proof of concept for a novel approach that quantitatively links the effects of altered molecular pathways onto neuronal population dynamics. As a novelty, we connect chemical compounds measured with positron emission tomography (PET) with neural function in TVB addressing the phenomenon of hyperexcitability in AD related to the protein amyloid beta (Abeta). We construct personalized virtual brains based on an averaged healthy connectome and individual PET derived distributions of Abeta in patients with mild cognitive impairment (MCI, N = 8) and Alzheimer's Disease (AD, N = 10) and in age-matched healthy controls (HC, N = 15) using data from ADNI-3 data base (http://adni.loni.usc.edu). In the personalized virtual brains, individual Abeta burden modulates regional Excitation-Inhibition balance, leading to local hyperexcitation with high Abeta loads. We analyze simulated regional neural activity and electroencephalograms (EEG).

Results: Known empirical alterations of EEG in patients with AD compared to HCs were reproduced by simulations. The virtual AD group showed slower frequencies in simulated local field potentials and EEG compared to MCI and HC groups. The heterogeneity of the Abeta load is crucial for the virtual EEG slowing which is absent for control models with homogeneous Abeta distributions. Slowing phenomena primarily affect the network hubs, independent of the spatial distribution of Abeta. Modeling the N-methyl-D-aspartate (NMDA) receptor antagonism of memantine in local population models, reveals potential functional reversibility of the observed large-scale alterations (reflected by EEG slowing) in virtual AD brains.

Discussion: We demonstrate how TVB enables the simulation of systems effects caused by pathogenetic molecular candidate mechanisms in human virtual brains.

Introduction

Neurodegenerative diseases (NDD) gain increasing socioeconomic relevance due to an aging society (WHO, 2011; Wimo et al., 2011, 2017; Xu et al., 2017). The Alzheimer's Association's latest report estimates the yearly cost of Alzheimer's disease (AD) treatment in the U.S. at $277 billion (Alzheimer's Association, 2018). By 2050 this number is expected to rise as high as $1.1 trillion. According to the report, early diagnosis could save up to $7.9 trillion in cumulated medical and care costs by the year 2050. While the prevalence of AD—the most common cause of dementia and the most common NDD in general—increases, its cause is still not understood, nor is there a cure. Our understanding of their pathogenesis and classification remain insufficient. Therefore, we aim to integrate clinical data from molecular biology and neurology, using nonlinear systems theory. Our aim is to build predictive models for health-outcome and cognitive function by individual virtual brain simulations using The Virtual Brain (TVB; thevirtualbrain.org) platform (Ritter et al., 2013; Sanz Leon et al., 2013). TVB integrates various empirical data in computational models of the brain that allow for the identification of neurobiological processes that are more directly linked to the causal disease mechanisms than the measured empirical data. Biomedical sciences are currently lacking a mapping between the degree and facets of cognitive impairments, biomarkers from high-throughput technologies, and the underlying causal origins of NDD like AD. The imperative for the field is to identify the features of brain network function in NDD that predict whether a person will develop dementia. The heterogeneity of NDD makes it difficult to develop robust predictions of cognitive decline. This can be addressed by large prospective studies where there is potential for participants to develop NDD. It is difficult in general to predict individual disease progression and this is a particular challenge in complex nonlinear systems, like the brain, where emergent features at one level of organization (e.g., cognitive function) can come about through the complex interaction of subordinate features (e.g., network dynamics, molecular pathways, gene expression). The Virtual Brain takes into account the principles of complex adaptive systems and hence poses a promising tool for identifying mechanistic predictive biomarkers for NDD. Due to the high dimensionality of brain models and the even greater complexity of the to-be-simulated brain states, selecting the used modeling approach carefully for a specific question of interest is essential.

The candidate biological mechanism under investigation in the present study is related to amyloid beta (Abeta), a protein that is an oligomeric cleavage product of the physiological amyloid precursor protein (APP) (Bloom, 2014; Selkoe and Hardy, 2016). The soluble oligomers have the tendency for polymerization (Sadigh-Eteghad et al., 2015; Selkoe and Hardy, 2016). Due to their non-physiological configuration they aggregate and accumulates in brain tissue—a process that starts already in early preclinical stages of AD, i.e., many years before the onset of symptoms—typically in the fifth decade of life (Braak and Braak, 1997)—as shown in rodent models (Busche et al., 2012) and human studies (Klunk et al., 2007; Jack et al., 2009). Aggregated Abeta and its intermediates, soluble Abeta oligomers, can act directly neurotoxic (Hardy and Selkoe, 2002; Prasansuklab and Tencomnao, 2013; Selkoe and Hardy, 2016) and have been found intra- or extra-cellularly (Hardy and Selkoe, 2002; Walsh and Selkoe, 2007; Selkoe and Hardy, 2016). Those findings led to the hypothesis that the deposition of Abeta poses an initial step in the pathology of AD while Abeta has been suggested as a key feature in the pathogenesis of AD leading to major changes in the functionality and structure of the brain (Klunk et al., 2007; Jack et al., 2009; Villemagne et al., 2009). The goal of the present study is to incorporate the hypothesized qualitative and quantitative effects of Abeta on neuronal population dynamics into our brain network models, i.e., adding mathematical models that describe how molecular changes alter population activity—so called cause-and-effect models. We will focus here on the disrupted inhibitory function of interneurons and consecutive hyperexcitability caused by Abeta—while we are aware of various other factors with potential roles for AD etiology, such as vascular changes (Love and Miners, 2016; Storck and Pietrzik Claus, 2018; Bannai et al., 2019), neuroinflammation (Heneka et al., 2015a,b; Wang and Colonna, 2019; Zhou et al., 2019), genetics (Mahley, 2016; Hudry et al., 2019; Takatori et al., 2019), environmental factors (Alonso et al., 2018; McLachlan et al., 2019) and concomitant proteinopathies others than Abeta pathology (Robinson et al., 2018a,b). Beside Abeta there is a second molecular hallmark associated with the pathogenesis of AD: the phosphorylated Tau “tubulin-associated unit” protein (Bloom, 2014; Guo et al., 2017; Tapia-Rojas et al., 2019) which contributes to microtubule stability in the neural cytoskeleton (Guo et al., 2017). One major argument in favor of the more prominent involvement of Abeta in the pathogenesis of AD, in contrast to Tau, is its higher specificity to AD and its appearance in the early familial variants of AD, where the molecular pathway is better understood (Blennow et al., 2006; Klunk et al., 2007; Villemagne et al., 2009). Therefore, most therapeutic strategies in the past targeted Abeta. Yet recently three clinical trials with antibodies against Abeta had to be terminated in phase III: aducanumab (Biogen, 2019; Chiao et al., 2019), crenezumab (Salloway et al., 2018; Roche, 2019), and solanezumab (Doggrell, 2018; Honig et al., 2018) did not meet the expectations to act in a disease-modifying manner slowing down the cognitive decline (Selkoe and Hardy, 2016). Nevertheless, there are still studies ongoing, e.g., with BAN-2401, an antibody against soluble monomeric Abeta protofribrils (as aducanumab) (Logovinsky et al., 2016; Osswald, 2018; Panza et al., 2019). A relevant percentage of clinically diagnosed AD patients show additional brain pathologies beside Abeta and Tau in autopsy (Robinson et al., 2018a). Even in the cases of neuropathological AD diagnosis (i.e., secured Abeta and Tau pathology in histology), 55% of cases also exhibited a pathology of alpha synuclein (which we would expect in synucleinopathies like Parkinson's disease) and 40% showed transactive response DNA binding protein 43kDa (TDP-43), a protein which we would expect in frontotemporal dementia or amyotrophic lateral sclerosis (Robinson et al., 2018b). Brain tissue of people who did not had relevant neurodegenerative brain changes in histological exams after death were showing Abeta in 50% and Tau pathology in 93% of the cases when using sensitive immunohistochemistry methods (Robinson et al., 2018b). Although Abeta and Tau are widely accepted as involved parts in the pathogenesis of AD and also define the disease entity (Jack et al., 2018), it remains unclear if they might be only epiphenomena of other contributing factors. This study hypothesizes a mechanistic role of Abeta in the disease process and builds a link between the molecular pathway alteration that leads to Abeta phenomenon of disinhibition and neural slowing in EEG (Figure 1). Our mechanistic modeling approach can help to understand the complex inter-dependencies between the involved factors in AD and will improve through iterative refinement.

FIGURE 1

Figure 1. Biology-infered cause-and-effect model: alteration of the molecular Abeta pathway in AD cause hyperexcitation in the neural mass model. An altered pathway from soluble Abeta monomers to oligomers to insoluble plaques leads to potentially neurotoxic Abeta accumulation (Hardy and Selkoe, 2002; Prasansuklab and Tencomnao, 2013; Selkoe and Hardy, 2016) that can be quantified by PET. Region specific Abeta burden leads to disinhibition in the neural mass model (Busche et al., 2008; Hazra et al., 2013; Ripoli et al., 2014; Ulrich, 2015; Ren et al., 2018)—thus building a bridge between molecular pathways and brain network modeling. For the evaluation of the used mathematical model, see the discussion section and Figures 12, 13. Parts of the figure are modified from Deco et al. (2017).

Near Abeta plaques, a shift in neural activity has been observed (Busche et al., 2008). In AD mouse models with overexpression of APP and Presenilin-1, the number of hyperactive neurons was increased near Abeta plaques. This shift in the neuronal activity was associated with decreased performances in memory tests. Neuronal hyperactivity could be reduced by GABA agonists, suggesting pathology due to impaired inhibition. In neocortical and dentate gyri, pyramidal cells have been found to increase network excitability in vivo in an AD mouse model with overexpression of Abeta, that led to membrane depolarization and increased firing rates. A study by Hazra et al. (2013) investigated an AD mouse model by stimulation of the perforant pathway. AD mice showed increased amplitude and larger spatial distribution of response after stimulation. The reason for this increased network excitability was due to impaired inhibitory neuron function, i.e., the inhibitory neurons of the molecular layer of the dentate gyrus in hippocampus were in part unable to produce action potentials, which resulted in a slower postsynaptic firing rate. Ulrich (2015) added Abeta to layer V pyramidal cells of rats. In their experiments they could show a decline in inhibitory postsynaptic currents (IPSCs), attributed to postsynaptic GABA_A receptor endocytosis after Abeta application. In a recent study by Ren et al. (2018) Abeta was found to increase excitability of pyramidal cells in the anterior cingulate cortex of mouse brain. The reason for hyperexcitability was again due to disturbed inhibitory input. Abeta seems to interact with the dopaminergic D1 receptor system. The D1 receptor regulates GABA release in fast-spiking (FS) inhibitory interneurons. By adding a D1 receptor antagonist to the cells they could reverse the effect of Abeta, increase IPSCs and decrease pyramidal excitability whereas D1 agonists had similar effects as Abeta. The underlying working model is that Abeta leads to dopamine release in dopaminergic neurons that activates D1 receptors at FS inhibitory interneurons and thus inhibits GABA release. As a consequence, the amplitude, frequency and total number of IPSPs is decreased. The instantaneous decrement of postsynaptic amplitude and frequency is also known as a toxic effect of Abeta in the glutamatergic system (Ripoli et al., 2014). Hence for the present modeling approach we decided to implement this Abeta dependent impaired inhibitory function. From the literature above, potential models for this disinhibiton could be either a lower IPSP amplitude or a lower firing rate or a combination thereof.

One already established drug that assesses the pathology of hyperexcitation is memantine, an N-methyl-D-aspartate (NMDA) antagonist. Memantine is recommended for the symptomatic treatment of severe AD as a mono- and combination therapy with cholinesterase inhibitors and should be also considered as possible treatment in moderate AD in the current version of the UK National Institute for Health and Care Excellence (NICE) guidelines of dementia management (Pink et al., 2018). However, normally it is considered as an alternative or addition to cholinesterase inhibitors (Pink et al., 2018). In contrast, memantine has shown in a current meta-analysis its efficacy to improve cognitive function and reduce behavioral disturbances in AD patients compared to placebo (Kishi et al., 2017). The effect was particularly caused by the moderate-to-severe AD patients (Chen et al., 2017; Kishi et al., 2017) and was also observable in combination therapies with acetyl cholinesterase inhibitors, with a significant superiority for the combination of memantine and donepezil compared to any cholinesterase monotherapies (Kishi et al., 2017). It therefore is also recommended as possible first-line therapy in AD (Kishi et al., 2017). In our study, we will evaluate “virtual memantine” interacting with the Abeta-derived hyperexcitation.

Changes in electroencephalography (EEG) are described in AD as a general and progressive slowing of brain oscillations. In AD, cognitive decline and ¹⁸F-fluorodeoxyglucose (FDG) PET signal decreases are linked with increased left temporal power in the delta and the theta frequency bands, whereas temporo-parieto-occipital alpha band coherence decreases and delta coherence increases (Loewenstein et al., 1989; Rice et al., 1990; Malek et al., 2017). Moreover, the spatial appearance of slow rhythms and hypometabolism in FDG PET have been linked (Dierks et al., 2000; Babiloni et al., 2016). A recent study produced similar findings in magnetoencephalography (MEG): A global increase of theta and a frontal increase of delta were correlated with entorhinal atrophy and glucose hypometabolism (Nakamura et al., 2018). In summary, a global slowing has been reported for AD, in particular a shift from alpha to theta and delta activity (Loewenstein et al., 1989; Rice et al., 1990; Dierks et al., 2000; Babiloni et al., 2016; Malek et al., 2017; Nakamura et al., 2018).

As a consequence of these findings, we will focus in our modeling approach on three main aspects of AD:

1. Spatial heterogeneous Abeta distribution in the brain

2. Hyperexcitation caused by impaired inhibitory function

3. Slowing of neural frequencies.

For Abeta, we propose a change in local neuronal excitability. Therefore, we construct a model of a healthy “standard brain” with an averaged structural connectivity (SC) with inferred micro-scale characteristics of excitation in those areas where a deposition of Abeta is found. We will infer this information about the local distribution of Abeta from individual AV-45 (florbetapir) positron emission tomography (PET) images. AV-45 is a PET tracer which binds to Abeta (Clark et al., 2011; Ossenkoppele et al., 2015; Morris et al., 2016; Martinez et al., 2017). Ante-mortem Abeta PET imaging can be related to post mortem Abeta pathology in brain tissue (Murray et al., 2015), corresponding to the THAL phases of Amyloid deposition (Thal et al., 2002)—as well as Tau PET (Schöll et al., 2016) can be related to the BRAAK stages of neurofibrillary tangles (Braak and Braak, 1991, 1997; Braak et al., 2006). This has led to updated diagnostic criteria for AD, wherein Abeta and Tau PET can be used equivalently to neuropathology for AD diagnosis (Jack et al., 2018).

We investigate three clinical diagnostic groups of age- and gender-matched healthy controls (HC), individuals with mild cognitive impairment (MCI) and AD patients [see method section Alzheimer's Disease Neuroimaging Initiative (ADNI) Database and Table 1]. For the simulated EEG and the underlying local neural activity frequency we expect a slowing in rhythms and particular a shift from alpha to theta activity with disease progression. Finally, we will simulate the effect of an anti-excitotoxic drug, the NMDA antagonist memantine for which we expect a reversal of the observed EEG slowing.

TABLE 1

Table 1. Basic epidemiological information of the study population.

We will in the following provide an overview of the fundamentals of the here employed brain simulation technique. The particular strength of computational connectomics (Ritter et al., 2013; Kringelbach et al., 2015; Deco et al., 2017) or brain network modeling (BNM) is to unite various kinds of information in a single biophysically plausible framework (Breakspear, 2017). BNM are typically structurally informed (or constrained) by (a) geometric information of the brain, e.g., via T1 magnetic resonance imaging (MRI), and (b) the structural connectivity (SC) derived from the tractography of diffusion MRI that is supposed to represent the white matter fiber tracts (Jirsa et al., 2002; Schirner et al., 2015). The static three-dimensional scaffold of the brain is brought to life through the implementation of mathematical models, which generate activity at each brain region or node of the network, the so-called neural masses or population models (Spiegler and Jirsa, 2013; Sanz-Leon et al., 2015; Cabral et al., 2017). Population models are reduced descriptions of microscopically detailed neuronal networks (Wilson and Cowan, 1972; Zetterberg et al., 1978; Hindmarsh and Rose, 1984; Jansen and Rit, 1995; Wong and Wang, 2006; Stefanescu and Jirsa, 2008; Sanz-Leon et al., 2015)—inferred for example with methods of mean field theory (Deco et al., 2008; Jirsa, 2009; Bojak et al., 2010). They describe the so called meso-scale of the brain (Deco et al., 2008; Wright and Liley, 2010), i.e., population activity as captured with imaging methods like EEG, MEG and fMRI. Some neural mass models (NMM) are linked to (and still reflect to a certain degree) neurophysiological processes at the microscopic scale while others mathematically describe the observed lumped biological behavior not differentiating between underlying neurophysiological processes (phenomenological models). Time delays in the interaction between nodes (Jirsa and Kelso, 2000; Jirsa et al., 2002; Spiegler and Jirsa, 2013; Sanz-Leon et al., 2015) are critical for the spatiotemporal organization of the evolving activity patterns in the brain (Petkoski et al., 2016, 2018). Measured functional brain data such as EEG, MEG or functional MRI (fMRI) are used to tune the mathematical models—i.e., to fit selected free parameters of the model—to faithfully reproduce selected empirical features (Honey et al., 2007; Ghosh et al., 2008; Sotero and Trujillo-Barreto, 2008; Bojak et al., 2010; Jirsa et al., 2010; Ritter et al., 2013; Sanz-Leon et al., 2015; Kunze et al., 2016). By performing a systematic model parameter exploration, using e.g., brute force exhaustive parameter space searches, Monte-Carlo methods or weighted optimization algorithms, we can identify the optimal parameter configuration to portray the empirical functional phenomena. Thereby, we obtain indices of the brains individual function in relation to the explored parameters. This approach opens various possibilities to not only describe dependencies (i.e., correlations), but to make statements about potential underlying causal processes, i.e., mechanisms.

In this study we used TVB, an open source neuroinformatics platform (Ritter et al., 2013; Sanz Leon et al., 2013; Sanz-Leon et al., 2015; Stefanovski et al., 2016) (www.thevirtualbrain.org) for large-scale BNM simulations. We have already established the software TVB, and applied it to normative datasets, stroke, epilepsy, brain tumors, and neurodegenerative disease. For example, in stroke recovery, TVB models of patients were built using the patient's structural neuroimaging data, and the dynamics of local populations were tuned to fit the patient's functional neuroimaging data (Falcon et al., 2015, 2016). The obtained parameters for excitatory/inhibitory (EI) balance of local neuronal populations predicted the patient's response to rehabilitation up to 1 year after therapy. Our work on epilepsy was able to infer seizure propagation with a model based on the patient's own diffusion weighted MRI and stereotaxic EEG (Jirsa et al., 2017; Proix et al., 2017). Moreover, positive surgical outcome was strongly associated with the epileptogenic zone that was excised as predicted by the patient's TVB model. Previous work with AD patients (n = 16), controls (n = 73), and persons with amnestic MCI (n = 35), all from the Sydney Memory and Aging Study, confirms the benefit of using the model parameters to characterize cognitive status (Zimmermann et al., 2018).

TVB provides several types of NMMs. In the present study, we selected a NMM that can simulate EEG and enables us to implement disinhibition. The wiring pattern of cortical circuitry is characterized by recurrent excitatory and inhibitory loops, and by bidirectional sparse excitatory connections at the large-scale (Schüz and Braitenberg, 2002). Several NMMs therefore feature projection neurons aka pyramidal cells with long axons projecting to distant cortical regions and local excitatory and inhibitory feedbacks (Lopes da Silva et al., 1974; Freeman, 1975; Jansen and Rit, 1995). The NMM by Jansen-Rit comprises an elementary circuit of three interconnected NMMs (Figure 2) describing a cortical area (or column). It has been used to explain both epilepsy-like brain activity (Wendling et al., 2000, 2002) and various narrow band oscillations ranging from the delta to the gamma frequency bands (David and Friston, 2003) including intracranial EEG (Spiegler and Jirsa, 2013). The Jansen-Rit model has been explored extensively on a single population level (Wendling et al., 2002; David and Friston, 2003; Spiegler et al., 2011) and in BNMs (Sotero et al., 2007; Merlet et al., 2013; Kunze et al., 2016). The Jansen-Rit model has a rich dynamic repertoire, which was extensively described before (Spiegler et al., 2010).

FIGURE 2

Figure 2. Postulated Abeta effect and its implementation to the Jansen-Rit model. (A) The virtual brains are based on averaged healthy connectomes and constrained by the individual regional burden of Abeta [figure modified from Deco et al. (2017)]. (B) In our simulation, increased excitability is caused by a shift in E/I balance, i.e., a slowed filter action in the transmission from inhibitory interneurons to pyramidal cells. In the background a histological representation of the cortical layers: excitatory pyramidal cells (υ₃) and excitatory interneurons (υ₁) are (exemplarily) located in layer V (internal pyramidal layer), while the inhibitory stellate (inter-neurons (υ₂) are located in layer IV (internal granular layer). In layer I (molecular layer) we see the dendrites of the pyramidal cells, where the input from the interneurons happens. The effect to the other neuron populations is represented by m₁₋₃ [background is a modified version of figure 13 from Schmolesky (2005)], license: https://creativecommons.org/licenses/by-nc/4.0/). (C) Schematic illustration of the three interacting neural masses in the Jansen-Rit population model. The reduced inhibition is mediated by negative influence of the local Abeta burden on the inhibitory time constant τ_i (see main text for more detailed explanation). This is intended to lead to an increased activity and higher output of the pyramidal cell population. The excitatory impulse response function (IRF) is specified as h_e(t) = tH_e exp(–t/τ_e)/τ_e, the inhibitory IRF is specified as h_i(t, β) = tH_i exp(–t/τ_i(β))/τ_i(β) (Equations 1, 2). These IRFs can be translated into second-order ordinary differential equations, see Equations 3–5. For explanation of the used variables, see Table 2 [figure modified from Spiegler et al. (2010)]. (D) Virtual EEG as the simulation output (projection of oscillating membrane potentials to the scalp surface) reveals a shift from alpha to theta activity in AD participants. Shown is a 5 second period of exemplary EEG channel at location T7 in participant 21 (HC, above) and 4 (AD, below). The ordinate is showing the dimensionless correlate for electric potential Φ. The exemplary timeseries shows a typical simulation result in the study: in the alpha mode, which was the starting point of the Jansen-Rit model without the effect of Abeta, it produces monomorphic alpha activity with amplitude modulations (above). Mainly exclusively in the AD virtual brains a much more irregular theta rhythm appears (below).

TABLE 2

Table 2. Used parameters for each Jansen-Rit element in the large-scale brain network (Jansen and Rit, 1995).

Specifically we chose the Jansen-Rit model for the present study due to the following considerations:

(1) The Jansen-Rit model comprises three interacting neural masses (representing different cellular populations) in each local circuitry: pyramidal cells, inhibitory, and excitatory interneurons (Figure 2B). This is unique and opens the possibility to simultaneously model disinhibition, i.e., an impairment of the inhibitory neuronal subpopulation in one neural mass, and an anti-NMDAergic effect, i.e., a downscaled transmission from excitatory interneurons to pyramidal cells, at the same time.

(2) The ratio of excitatory and inhibitory time constants τ_e/τ_i in the Jansen-Rit model is suitable to model the effect of Abeta on the inhibitory interneurons (by affecting the transmission from inhibitory interneurons to pyramidal cells, Figures 2B,C) and is also known to have an effect on the simulated neural frequency (Wendling et al., 2002; Spiegler et al., 2010). Because oscillations emerge in the Jansen-Rit model of a brain region due to the interplay of positive and negative feedback loops (excitatory and inhibitory interneurons), a change in one of the time constants does not necessarily slow down or speed up rhythms. However, if both excitatory and inhibitory time constants, τ_e and τ_i are scaled simultaneously and uniformly, the local equilibrium of interaction between the neural masses remains the same but the time signature such as frequency changes [see Figure 9 in Spiegler et al. (2010), and see Chapter 3.1.3 “Model Equivalence” and Chapter 3.1.4 “Normalization” in Spiegler (2011)]. To conclude, higher τ_i does not necessarily lead to slower rhythms and vice versa.

(3) Jansen-Rit can simulate physiological rhythms observable in local field potentials (intracranially), stereo-EEG (sEEG), scalp EEG, and MEG (Jansen and Rit, 1995; Spiegler et al., 2010; Sanz-Leon et al., 2015).

Our hypothesized effect of local Abeta deposition as inferred from subject-specific AV-45 PET is a decrease of local inhibition (Busche et al., 2008; Grienberger et al., 2012; Limon et al., 2012; Verret et al., 2012; Hazra et al., 2013; Ripoli et al., 2014; Ren et al., 2018), which leads to a relatively stronger local excitation. This theory allows us translation of the Abeta distribution into the altered dynamics of a population model (Equation 14 and Figure 3). We use an averaged healthy SC to control the effect of individual differences in the connectome. I.e., in our simulations the distribution of Abeta is the only individual factor and can therefore be seen as the cause of any differences between the participants. The hypothesized microscale (synaptic), spatially distributed effect is assumed to develop an effect at the population (mesoscale) level and to eventually propagate to the large-scale of the whole brain. A schematic illustration of this concept is provided in Figures 1, 2.

FIGURE 3

Figure 3. Graphs of the sigmoid transfer function of Abeta. The abscissa represents the Abeta burden β_a, the ordinate represents the reciprocal S_τ(β_a) of the inhibitory time constant τ _i. See Equation 14.

Methods

Alzheimer's Disease Neuroimaging Initiative (ADNI) Database

Empirical data 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. The primary goal of ADNI has been to test whether serial MRI, PET, other biological markers, and clinical and neuropsychological assessment can be combined to measure the progression of MCI and early AD. For up-to-date information, see http://www.adni-info.org.

In the presently ongoing trial, ADNI-3, the measurements contain T1, T2, DTI, fMRI, Tau PET, Abeta PET, and FDG PET for the participants. The total population of ADNI-3 will contain data of about 2,000 participants (comprising AD, MCI, and HC, see http://adni.loni.usc.edu/adni-3/). As inclusion criterion for AD patients the diagnosis criteria of NINCDS-ADRDA from 1984 were used, which contains only clinical features (McKhann et al., 1984). Inclusion criteria for both HC and MCI were a Mini Mental State Examination (MMSE) score between 24 and 30 as well as age between 55 and 90 years. For MCI in addition, the participant must have a subjective memory complaint and abnormal results in another neuropsychological memory test. To fulfill the criteria for AD, the MMSE score had to be below 24 and the NINCDS-ADRDA criteria for probable AD had to be fulfilled (McKhann et al., 1984). Imaging and biomarkers were not used for the diagnosis. For the full inclusion criteria of ADNI-3 see the study protocol (page 11f in http://adni.loni.usc.edu/wp-content/themes/freshnews-dev-v2/documents/clinical/ADNI3_Protocol.pdf). An overview of the epidemiological characteristics of the participants included in this study can be found in Table 1.

Data Acquisition and Processing

All images used in this study were taken from ADNI-3. To reach comparable datasets, we used only data from Siemens scanners with a magnetic field strength of 3T (models: TrioTim, Prisma, Skyra, Verio). However, some acquisition parameters differed slightly. See Supplementary Material with Supplementary Tables 1–6 for the metadata. The following imaging modalities were included: T1 MPRAGE. TE = 2.95–2.98 ms, TR = 2.3 s, matrix and voxel size differ slightly. FLAIR. TE differs slightly, TR = 4.8 s, matrix size = 160 · 256 · 256, voxel size differs slightly. DWI (only for 15 HC participants to create an average healthy SC). TE = 56–71 ms, TR = 3.4–7.2 s, matrix size = 116 · 116 · 80, voxel size = 2 · 2 · 2, bvals = [0, 1000] or [0, 500, 1000, 2000], bvecs = 49 or 115. Siemens Fieldmaps and PET Data (AV-45 for Abeta). The preprocessing of imaging data can be subdivided in that of structural images, DWI, and PET.

Structural MRI

We calculated an individual brain parcellation for each included participant of ADNI-3. We followed the minimal preprocessing pipeline (Glasser et al., 2013) of the Human Connectome Project (HCP) for our structural data using Freesurfer (Reuter et al., 2012) (https://surfer.nmr.mgh.harvard.edu/fswiki/FreeSurferMethodsCitation), FSL (Smith et al., 2004; Woolrich et al., 2009; Jenkinson et al., 2012) and connectome workbench (https://www.humanconnectome.org/software/connectome-workbench). Therefore, we used T1 MPRAGE, FLAIR and fieldmaps for the anatomical parcellation and DWI for tractography. This consists of a Prefreesurfer, Freesurfer, and Postfreesurfer part. We skipped the step of gradient non-linearity correction, since images provided by ADNI already are corrected for this artifact. Also, the MNI templates were used at 1 mm resolution instead of 0.7 mm. In the Freesurfer pipeline we skipped the step of downsampling our data from 0.7 to 1 mm³, and all recon-all and intermediate steps were performed with the original image resolution. We then registered the subject cortical surfaces (32 000 vertices) to the cortical parcellation of Glasser et al. (2016) using the multimodal surface matching (MSM, see Robinson et al., 2014) tool. For the registration we used cortical thickness, MyelinMaps, cortical curvature and sulc from the subject and template surface. We mapped the parcellation on the surface back into the gray matter volume with connectome workbench. This volume parcellation surfed as the mask for the connectome and PET intensity extraction.

PET Images

We used the preprocessed version of AV-45 PET. These images had following preprocessing already performed by ADNI: Images acquired 30–50 min post tracer injections: four 5-min frames (i.e., 30–35 min, 35–40 min…). These frames are co-registered to the first and then averaged. The averaged image was linearly aligned such that the anterior-posterior axis of the subject is parallel to the AC-PC line. This standard image has a resolution of 1.5 mm cubic voxels and matrix size of 160 · 160 · 96. Voxel intensities were normalized so that the average voxel intensity was 1. Finally, the images were smoothed using a scanner-specific filter function. The filter functions were determined in the certification process of ADNI from a PET phantom. We used the resulting image and applied the following steps: Rigid aligning the PET image to participants T1 image (after being processed in the HCP structural pipeline). The linear registration was done with FLIRT (FSL). The PET image was than masked with the subject specific brainmask derived from the structural preprocessing pipeline (HCP). To obtain the local burden of Abeta, we calculated the relative intensity to the cerebellum as a common method in the interpretation of AV-45-PET, because it is known that the cerebellum does not show relevant AV-45 PET signals and can therefore act as a reference region for inter-individual comparability between patients (Clark et al., 2011; Meyer et al., 2018). The intensity of gamma radiation, which is caused by a neutralization reaction between local electrons and the emitted positrons of the nuclear tracer is measured for each voxel in the PET image and divided to the cerebellar reference volume: the standardized uptake value ratio (SUVR). We therefore receive in each voxel a relative Abeta burden β which is aggregated according to the parcellation used for our present modeling approach (see below). Thus, we obtain a value β_a for the Abeta burden in each brain region a. The cerebellar white matter mask was taken from the Freesurfer segmentation (HCP structural preprocessing). The image was then partial volume corrected using the Müller-Gärtner method from the PETPVC toolbox (Thomas et al., 2016). For this step the gray (GM) and white matter segmentation from Freesurfer (HCP structural preprocessing) was used. Subcortical region PET loads were defined as the average SUVR in subcortical GM. Cortical GM PET intensities were mapped onto the individual cortical surfaces using connectome workbench tool with the pial and white matter surfaces as ribbon constraints. Using the multimodal parcellation from Glasser et al. (2016) we derived average regional PET loads.

DWI

We calculated individual tractography only for included HC participants of ADNI-3 to average them to a standard brain template (see section Virtual Human Standard Brain Template Out of Averaged Healthy Brains below). Preprocessing of the diffusion weighted images was mainly done with the programs and scripts provided by the MRtrix3 software package (http://www.mrtrix.org).

The following steps were performed:

Dwidenoise. Denoising the DWI data using the method described in Veraart et al. (2016).

Dwipreproc. Motion and eddy current correction using the dwipreproc wrapper script for FSL (https://mrtrix.readthedocs.io/en/latest/dwi_preprocessing/dwipreproc.html).

Dwibiascorrect. B1 field inhomogeneity correction using ANTS N4 algorithm.

Diw2mask. brainmask estimation from the DWI images.

Dwiintensitynorm. DWI intensity normalization for the group of participants.

Dwi2response. The normalized DWI image was used to generate a WM response function. We used the algorithm described by Tournier et al. (2013) in this step.

Average_response. An average response function was created from all participants.

Dwi2fod. Using the spherical deconvolution method described by Tournier et al. (2007) we estimated the fiber orientation distribution using the subject normalized DWI image and the average response function. From the DWI data a mean-b0 image was extracted and linear registered to the T1 image. The inverse of the transform was used to bring the T1 brain masked and aparc+aseg image (from HCP structural preprocessing) into DWI space. The transformed aparc+aseg image was used to generate a five tissue type image.

Tckgen. Anatomical constrained tractography (Smith et al., 2012) was performed using the iFOD2 algorithm (Tournier et al., 2010). Tracks in the resulting image were weighted using SIFT2 algorithm (Smith et al., 2015). We mapped the registered parcellation from Glasser back into the volume. The cortical and subcortical regions than were used to merge the tracks into a connectome.

EEG Forward Solution in TVB

After structural preprocessing with the HCP pipeline we used the individual cortical surfaces and T1 images to compute the person specific Boundary Element Models in Brainstorm (Tadel et al., 2011). Scalp, outer, and inner skull were modeled with 1922 vertices per layer. Using the default “BrainProducts EasyCap 65” EEG cap as locations for the signal space and the cortical surface vertices as source space. The leadfield matrix was estimated using the adjoint method in OpenMEEG with the default conductivities 1, 0.0125 and 1 for scalp, skull and brain, respectively. Because we are performing region-based simulations only (i.e., no vertex-wise modeling), the leadfield matrix was simplified by summing the coefficients of vertices that belong to the same region. EEG signal was generated by matrix multiplication of the neural time series with the lead field matrix.

Virtual Human Standard Brain Template Out of Averaged Healthy Brains

We use the SCs of all ADNI-3 participants of the group HC, derived from the diffusion-weighted and structural MRI, to average them to one connectome matrix. Two of the HC participants included in the average template were excluded for simulations because it was impossible to compute their leadfield matrices for EEG calculation. Therefore, we use an arithmetic mean C_μ = (${\sum }_i^n$ _{= 1} C_i)/n = (C₁ + C₂ + … + C_n)/n, wherein C_μ is the averaged SC matrix, n is the number of HC participants and C_i is the individual SC matrix. The SC matrix and the organization of the corresponding graph can be found in Figure 4. As it can be seen in Figure 4B, general characteristics of physiological SCs as symmetry, laterality, homology, and subcortical hubs are maintained in the averaged connectome. By choosing an averaged SC instead of individual SCs, it was possible to control all factors except of the individual Abeta distribution supporting our intention to compare the simulated activity that resulted from a “pathogenic” modification by Abeta.

FIGURE 4

Figure 4. Underlying average HC structural connectome. (A) SC Matrix of the underlying averaged SC, showing the DWI-derived connections weights. Three hundred and seventy nine regions are in following order: 180 left cortical regions, 180 right cortical regions of Glasser parcellation (Glasser et al., 2016), 9 left subcortical regions, 9 right subcortical regions, 1 brainstem region. It gets obvious the difference between interhemispherical commissural fibers (lower weights, with a slightly pronounced diagonal between homologous regions) and intrahemispherical association fibers (higher weights). Moreover, we can observe the strong connection pattern of the subcortical areas (above region 360). (B) Graph of the underlying SC. As a threshold, only the strongest 5% of connections were kept for binary transformation to the adjacency matrix for the graph. Node positions are derived from the inner structure of the graph by a “force” method (Fruchterman and Reingold, 1991), assuming stronger forces and therefore smaller distances between tightly connected nodes. It can be seen that the laterality is kept in the graph structure (also for subcortical regions) and the whole graph is highly symmetric. Node size linearly represents the graph theoretical measure of structural degree for each node. Most important hubs are subcortical regions. The shown features of symmetry, laterality, homology, and subcortical hubs indicate that the averaged SC still kept its physiological characteristics.

Cause-and-Effect Model of Abeta in the Jansen-Rit Model

The dynamics of the Jansen-Rit model show a rich parameter dependent behavior (Spiegler et al., 2010). A bifurcation analysis of the single population Jansen-Rit model (in contrast to network embedded interacting populations) catalogs and summarizes the repertoire of the model. Bifurcation here refers to a qualitative change in the system behavior with respect to parameter changes. Qualitative changes can be for instance the shift from waxing and waning alpha rhythm as observed in resting human brains to spike wave discharges as observed during epileptic seizures. Bifurcation diagrams explore the qualitatively different states [divided by bifurcations, see Supplementary Figure 1, from Spiegler et al. (2010)]. The bifurcation analysis revealed an important feature of the Jansen-Rit model, which is bistability, that is, the coexistence of two stable states for a certain parameter range (i.e., regime). The bistable regime allows the coexistence of two self-sustained oscillatory states for the standard parameter configuration (Jansen and Rit, 1995) and Table 2 of which one state generates rhythmic activity in the alpha band and the other one produces slower big spike-wave complexes in theta rhythm. Changes in the kinetics of excitatory and inhibitory PSPs (i.e., changes of time constants) change the model behavior in a way which makes it suitable to scale, that is, to speed up or to slow down local dynamics (Spiegler et al., 2010)—and therefore to scale the global frequency, too. The results of the systematic parameter exploration of the excitatory and inhibitory time constants are summarized in Supplementary Figure 2. For our study, we constructed the model to portrait a wide range of physiological neural rhythms by using a fast limit cycle with alpha and beta frequencies and a slow limit cycle with theta and delta frequencies. To achieve this dynamic behavior of two limit cycles, we used first a very low input on the pyramidal cells (firing rate 0.1085/ms) and no input on the inhibitory interneurons to not overlay the Abeta effects we introduce here. Therefore in the “healthy” condition without the effect of Abeta, the system operates near the subcritical Andronov-Hopf and the saddle-saddle bifurcations (leftmost region in Supplementary Figure 1). For the time constants, we used the area of alpha rhythm (blue area in Supplementary Figure 2) as control condition without any effect of Abeta. The detailed parameter settings can be found in Table 2.

The information about the local Abeta burden is derived from the individual AV-45 PET. As there exists no established clinical standard for SUVR cut-off thresholds differentiating normal form pathological Abeta loads. To scale the possible neurotoxic effect in a realistic way, we need to approximate at what point Abeta toxicity occurs. Following the literature, a 96% correlation to autopsy after Abeta PET was achieved via visual assessment of PET images. The corresponding SUVR cut-off was 1.2 (Clark et al., 2011). Another study showed a higher cut-off point at SUVR ≥ 1.4 for a 90% sensitivity of clinically diagnosed AD patients with an abnormal Abeta PET scan (Jack et al., 2014). We use here the higher cut-off threshold of SUVR ≥1.4. Consequently, we propose a cause-and-effect model for Abeta that is mapping molecular changes to computational brain network models:

The inhibitory time constant τ_i in each point is a function of β_a. The higher Abeta SUVR, the higher τ_i and therefore the filter action for the synaptic transmission is slower. We decided for this implementation via a synaptic filter slowing because of several reasons:

1. We are focusing on disease linked alterations of EEG frequencies. Hence, we intended to assess a model feature that is already known to be frequency-effective, i.e., it can vary resulting simulated EEG frequencies. From former explorations of the Jansen-Rit-model we know that the neural frequencies are influenced by the ratio of excitatory and inhibitory time constants (Spiegler et al., 2010).

2. Cellular studies are supporting the hypothesis of altered inhibition as a cause for hyperexcitation (Hazra et al., 2013; Ripoli et al., 2014; Ren et al., 2018)—hence we decide for mapping Abeta on the inhibitory time constant leading to a disturbed Excitatory-Inhibitory (E/I) balance.

3. By using a time-effective feature, we intended to differentiate the micro-scale neurotoxic effect of Abeta on synaptic level (Ripoli et al., 2014; Ulrich, 2015; Ren et al., 2018) from connectivity-effective phenomena on a larger scale, which could e.g., be modeled by an alteration of connection strength.

A detailed exploration of the effects that we introduce by this model can be found in the discussion section.

We develop a transform function to implement the PET SUVR in parameters of the brain network model. Specifically, we postulate a sigmoidal decrease function that modifies the default value for inhibitory time constant τ_i (Equation 14 and Figure 3). We assume the healthy brain without super-threshold Abeta burden operates in a region of the parameter space, which is close to a network criticality. A criticality describes an area in the parameter space, where subtle changes of one variable can have a critical impact on others (Strogatz, 2015) (in this case bifurcations, see Supplementary Figure 1. The thresholding “cut-off” value β_off–differentiating normal form pathological Abeta burden—was chosen according to the literature, stating that only after a certain level of tracer uptake a region is considered pathological (β_off = 1.4, see above). The maximum possible Abeta burden value β_max was chosen to be the 95% percentile of the Abeta regional SUVR distribution across all participants. The midpoint of the sigmoid was chosen such that it was half the way between β_off and β_max. The steepness was chosen such that the function converges to a linear function between β_off and β_max.

Brain Network Model Construction and Simulation

For the reasons stated in the above introduction, for our simulation approach we selected the Jansen-Rit model (Zetterberg et al., 1978; Jansen and Rit, 1995; Wendling et al., 2000; David and Friston, 2003; David et al., 2006; Spiegler et al., 2010, 2011; Sanz-Leon et al., 2015; Kunze et al., 2016). The differential equations are presented in Equations 3–5 (Jansen and Rit, 1995). The employed parameter values can be found in Table 2.

As a general approach, the impulse response function (IRF) of a neural mass allows to transform an incoming action potential into a PSP by using a linear time-invariant system. The IRF is the transfer function of the system, which is convoluted with the incoming input (action potentials) to calculate the output (PSPs). The general form of the IRF is the systems output to a (infinitesimal short and high) Dirac impulse and can be estimated experimentally by using short impulses or step functions (Lopes da Silva et al., 1974; Freeman, 1975).

The excitatory IRF h_e(t) is therefore specified as

$\begin{array}{ll}{h}_{\text{e}}(t)={\text{t}H}_{\text{e}}\text{exp}(-\text{t}/{\tau }_{\text{e}})/{\tau }_{\text{e}},& (1)\end{array}$

where τ_e is the excitatory time constant (the time until the PSP reaches its maximum), H_e is a coefficient of the PSP amplitude and t is time.

Similarly, the inhibitory IRF h_i(t, β) is specified as

$\begin{array}{ll}{h}_{\text{i}}(t,\beta )={tH}_{\text{i}}\text{exp}(-\text{t}/{\tau }_{\text{i}}(\beta ))/{\tau }_{\text{i}}(\beta ),& (2)\end{array}$

with the same variables as above. As we will describe below in detail, the inhibitory IR is a function of the spatially distributed Abeta burden β, which affects the time characteristics τ_i(β) and therefore h_i(t, β) is a function of time and space.

These IRFs can be translated into second-order ordinary differential equations by interpreting them as Green's functions. See Spiegler et al. (2010) for a detailed explanation of the dimensional reduction used here.

The differential equations that describe the network of three neural masses are now presented in Equations 3–5. The variables used for the simulations are listed in Table 2:

Excitatory projections υ₁ onto pyramidal cells at location a in discretized space (a = 1, 2, …, N: N = 379 regions):

$\begin{array}{ll}\{\begin{array}{l}\text{d}{\upsilon }_{1,a}(t)/\text{d}t={\upsilon }_{4,a}(t)\\ \text{d}{\upsilon }_{4,a}(t)/\text{d}t=H_{\text{e}}(m_{3\text{T},a}(t)+c_{\mathrm{31}}S(c_{\mathrm{13}}{\upsilon }_{3,a}(t)))/{\tau }_{\text{e}}\\ -2{\upsilon }_{4,a}(t)/{\tau }_{\text{e}}-{\upsilon }_{1,a}(t)/{\tau }_{\text{e}}^2\end{array}& (3)\end{array}$

Inhibitory projections υ₂ onto pyramidal cells at location a in space:

$\begin{array}{ll}\{\begin{array}{l}\text{d}{\upsilon }_{2,a}(t)/\text{d}t={\upsilon }_{5,a}(t)\\ \text{d}{\upsilon }_{5,a}(t)/\text{d}t=c_{32}{H}_{\text{i}}S(c_{\mathrm{23}}{\upsilon }_{3,a}(t))/{\tau }_{\text{i}}({\beta }_a)\\ -2{\upsilon }_{5,a}(t)/{\tau }_{\text{i}}({\beta }_a)-{\upsilon }_{2,a}(t)/{\tau }_{\text{i}}^2({\beta }_a)\end{array}& (4)\end{array}$

Projections υ₃ of pyramidal cells a in space:

$\begin{array}{ll}\{\begin{array}{l}\text{d}{\upsilon }_{3,a}(t)/\text{d}t={\upsilon }_{6,a}(t)\\ \text{d}{\upsilon }_{6,a}(t)/\text{d}t=S({\upsilon }_{30,a})(t)/{\tau }_{\text{e}}\\ -2{\upsilon }_{6,a}(t)/{\tau }_{\text{e}}-{\upsilon }_{3,a}(t)/{\tau }_{\text{e}}^2,\end{array}& (5)\end{array}$

wherein c₃₁, c₁₃, c₂₃ are the local connectivity weights between the three neural masses. Equation (4) shows the spatial dependency of the activity of inhibitory interneurons projected onto the pyramidal cells by τ_i(β_a).

Taking into account the biologically plausible configuration of the Jansen-Rit model shown in Figure 2, we use several mathematical simplifications to reduce the model's dimensionality without loss of generality. Taking into account that Equation (1) is valid for all excitatory input at dendrites irrespective of the source allows for using one single IRF at the pyramidal cells

$\begin{array}{ll}{h}_{\mathrm{31}}\equiv h_{3T}\equiv h_e,& (6)\end{array}$

and, thanks to linearity, translates the summation of excitatory postsynaptic potentials

$\begin{array}{ll}{\upsilon }_{31}+{\upsilon }_{3\text{T}}\equiv {\upsilon }_1& (7)\end{array}$

into a sum of incoming firing rate, that is, m_{3T, a}(t) + c₃₁ S(c₁₃ υ_{3, a}(t)) in Equation (3), describing the excitatory projections onto pyramidal cells υ₁. This simplification is without restrictions, simply exploits the linearity of the operators and reduces the dimensionality by 2. Furthermore, to adjust notation, the postsynaptic potentials caused by the inhibitory neural mass at pyramidal cells are denoted as

$\begin{array}{ll}{\upsilon }_{\mathrm{32}}\equiv {\upsilon }_2,& (8)\end{array}$

and its kernel is as

$\begin{array}{ll}{h}_{\mathrm{32}}\equiv h_i.& (9)\end{array}$

The projecting variable of one brain region at location a to other regions in the network is the sum of excitatory and inhibitory postsynaptic potentials at the local neural mass of pyramidal cells

$\begin{array}{ll}{\upsilon }_{30,a}(t)={\upsilon }_{1,a}(t)-{\upsilon }_{2,a}(t)& (10)\end{array}$

transferred into a firing rate using a sigmoid. The general form of this transfer function is

$\begin{array}{l}S(\lambda )=(S_{\text{max}}-S_{\text{min}})/(1+\text{exp}(r_{\lambda }({\lambda }_0-\lambda )))\\ +S_{\text{min}}:0<S_{\text{min}}<S_{\text{max}},& (11)\end{array}$

with, λ = υ, S_υ, _max = 2e₀ and S_υ, _min = 0 for the potential-to-firing-rate transfer.

Incoming mean firing rates m_{3T, a}(t) at the pyramidal cells at location a from other brain regions b = 1, 2, …, N, where N is the number of 379 regions are given by

$\begin{array}{ll}{m}_{3\text{T},a}(t)=m_{3\text{T},0}+G{\displaystyle \sum {}_b{w}_{a,b}S({\upsilon }_{30,b})(t),}& (12)\end{array}$

where m_{3T, 0} is baseline input m_{3T, 0} = const. for ∀t and all locations ∀a. The global coupling factor G is a coefficient of the incoming activity and therefore scales the connections w_{a, b} incoming at location a from all b provided by the SC. Because of this, G is the crucial factor that moderates the influence of the network to each neural mass and therefore mediates the difference between an uncoupled systems (G = 0) and a strong connected system.

In all populations, the state variable [υ₁, υ₂, υ₃]_a(t) are the mean membrane potentials and the derivatives thereof with respect to time t, namely [υ₄, υ₅, υ₆]_a(t) represent the mean currents.

To model how the local Abeta load β_a, measured by the Abeta PET SUVR is affecting the inhibitory time constant we introduce a transfer function (Figure 3). The primary assumption of this transfer function is a dependency of the E/I balance on the local Abeta concentration as described above. With higher concentration of Abeta, we assume dynamic changes in the inhibitory population of the neural mass that lead to local hyperexcitation. To model this inside the existing Jansen-Rit equations, potential candidate parameters would be H_i and τ_i as well as c₃₂. The coefficient H_i is not a suitable candidate because it has no direct physiological correlate. The coupling coefficient c₃₂ corresponds best to synaptic transmission from inhibitory to pyramidal cells and therefore can be mainly seen as a receptor surrogate. The time constant τ_i acts as a filter for IPSPs and correlates best with the evidence of decreased IPSP firing rate (Busche et al., 2008; Grienberger et al., 2012; Limon et al., 2012; Verret et al., 2012; Hazra et al., 2013; Ripoli et al., 2014; Ren et al., 2018) and is moreover well explored for the Jansen-Rit model (Wendling et al., 2002; Spiegler et al., 2010). As described above, we expect this transfer function to behave in an asymptotic way for Abeta concentrations below and above a specific range. We determined the lower border at β_{a, off} = 1.4 and the upper border at the 95th percentile in our data at β_{a, max} = 2.65. By exploring the effects of τ_i in a single region model, we determined the effective range 14.29 ms ≤ τ_i <50 ms. Based on this range we defined the following three-conditional linear function

$\begin{array}{ll}\{\begin{array}{ll}{\tau }_{\text{i}}={\tau }_{\text{i},\text{min}}& \text{for}{{\beta }_a\le {\beta }_a}_{,\text{off}}\\ {\tau }_{\text{i}}({\beta }_a)=m·{\beta }_a-c& \text{for}{{\beta }_a}_{,\text{off}}{<{\beta }_a<{\beta }_a}_{,\text{max}}\\ {\tau }_{\text{i}}={\tau }_{\text{i},\text{max}}& \text{for}{{\beta }_a\ge {\beta }_a}_{,\text{max}}\end{array}& (13)\end{array}$

wherein τ_i,min and τ_i,max are the maximum and minimum values for τ_i,

m = (τ_i,max-τ_i,min)/(β_a,max-β_a,off) = 28.6 and

c = m · β_a,off-τ_i,min = 25.7.

Since this function is not differentiable in β_a,off and β_a,max, we used the sigmoid function Equation (12) instead, which is continuous and differentiable. Moreover, a sigmoid can be interpreted as the cumulative (of a logistic distributed) activity acquired by the PET of a small brain volume (voxel) with a low spatial resolution of about 2.5 mm and above (Moses, 2011).

Therefore, the Abeta transfer function is defined as

$\begin{array}{ll}\{\begin{array}{l}{\tau }_{\text{i}}({\beta }_a)=S^{-1}({\beta }_a)\\ r_{\beta a}=2\text{ln}(S_{\text{max}}·1\text{s}-1)/({{\beta }_a}_{,\text{off}}-{{\beta }_a}_{,\text{max}})\\ {\beta }_0=({{\beta }_a}_{,\text{off}}{+{\beta }_a}_{,\text{max}})/2,\end{array}& (14)\end{array}$

wherein r_βa is the slope of the sigmoid, β₀ is the midpoint of the sigmoid and the coefficients are chosen to fit the conditions explained before. In this function, τ_i appears as its reciprocal value τ_i⁻¹ as it is implemented in the code of TVB. Because τ_i is a time in ms, the inverse of τ_i is a rate of potential change, and does not directly correspond to a firing rate. The Abeta load affects the inhibitory rate following a sigmoid curve. The rate ranges between S_min and S_max and the time constant ranges consequently between 1/S_max and 1/S_min.

To simulate the model using TVB, physical space and time are discretized. The system of difference equations is then solved using deterministic Heun's method with a time step of 5 ms. We used a deterministic method to avoid stochastic influences since the simulation was performed in the absence of noise.

The system was integrated for 2 min and the last minute was analyzed in order to diminish transient components in the time series due to the initialization and settle the system into a steady state.

We explore a range of 0 ≤ G ≤ 600 which provides an overview about the possible population level behaviors at different states of network coupling. Because the coupling factor G has a crucial influence on the external input on the neuronal populations, this allows different regions to operate in different dynamical regimes, as it can be seen in the bifurcation diagrams of Supplementary Figure 1. Global coupling factor G that was sampled between G = 0 (i.e., isolated regions) and G = 600 with a step size of ΔG = 3. The initial values were taken from 4,000 random time points for each state variable in each region. The length of 2 min for the simulations was chosen with the aim to diminish possible transient components due to the initialization of state variables at t = 0. For analysis we used only the second minute of the simulated signals. No time delays are implemented in the large-scale network interactions since they are not required for the emergence of the here evaluated features and setting them to zero increases reduces required computation resources.

Spectral Properties of the Simulated EEG

In TVB, we simulate EEG as a projection of the oscillating membrane potentials inside the brain via its electromagnetic fields to the skin surface of the head (Sanz-Leon et al., 2015) using the individual lead field matrices which take into account the different impedances of white matter, gray matter, external liquor space, pia and dura mater, the skull and the skin. Our lead-field matrices considered the impedances of three compartment borders: brain-skull, skull-scalp and scalp-air (Jirsa et al., 2002; Bojak et al., 2010; Litvak et al., 2011; Ritter et al., 2013). The postsynaptic membrane potential (PSP) considered for the projection is the one of the pyramidal cells, as they contribute the mayor part to potential changes in EEG (Kirschstein and Köhling, 2009). The PSP is calculated by summing the synaptic input from excitatory and inhibitory subpopulations to the pyramidal cells. The baseline PSP was derived as the mean PSP across time for every region. For the LFP or EEG peak frequency, we computed the power spectrum using the “periodogram” function of the Scipy python toolbox (Jones et al., 2001). From the spectrogram the “dominant rhythm” was identified as the frequency with the highest power.

Results

Abeta-Inferred Dynamics Lead to Individual Spectral Patterns

We analyzed the dominant frequency in the simulated EEG and regional neural signal (referred to as local field potential (LFP) (Figures 5G–J).

FIGURE 5

Figure 5. Spectral behavior in individuals of the different groups. (A–F) Selected timeseries and spectrograms. On the left: power spectral density of neural activity for an exemplary region with different values for G. Abscissa is frequency, ordinate is an estimate of the spectral power (dimensionless equivalent of amplitude per 1 Hz). Colors are representing the EEG frequency bands from delta to beta, indicated with Greek letters (note that this is regional neural activity, not EEG). Corresponding time series on the right: neural activity at single regions, each showing 5 s. Abscissa is time, ordinate is a dimension-less equivalent of the electric potential. (A) Shows an irregular, amplitude modulated alpha to beta rhythm, (B) an irregular theta with some delta and alpha inside. In (C) we can observe a monomorphic spike signal with a theta/delta frequency of 3 Hz and higher order harmonies. (D) Shows a monomorphic (high) alpha rhythm, (E) shows the zero-line with a continuous power spectrum. (F) Time series of 10 regions in a G area of hypersynchrony. We can see here the synchronized signals in theta rhythm and multiple harmonies of higher order in the spectrogram. (G–J) Four exemplary participants with different types of frequency behaviors along the range of coupling. Shown are the regional simulated dominant frequencies (y) along global coupling G (x) for individual exemplary participants 4, 8, 12, and 23. See Supplementary Table 7 for participant IDs. Color indicates the density of regions with the same coordinates. The sources of the timeseries on the left (A–F) are marked in the plots. (G) Irregular or chaotic rhythm with two clusters in alpha and theta. AD participant 8. (H) Chaotic behavior for lower G, then harmonic and hypersynchronization. AD participant 4. (I) Early zero-line, with monomorphic alpha activity at very low G. HC participant 23. (J) Harmonic to zero-line rhythm, with a G area of hypersynchrony in alpha and theta, depending on G.

We observed a physiologically looking irregular behavior with two frequency clusters in the alpha and in the theta spectrum (Figure 5G). This behavior is expressed in the area of lower global coupling G for all 10 AD participants and in 3 out of 8 MCI and 4 out of 15 HC participants. The irregular time series and the broad continuous frequency spectra (Figure 5B) of network regime in 0 < G < 150 are indicative for deterministic chaos. Such chaotic network regimes in a BNM have already been reported using the same local dynamic model [Figure 2 in Kunze et al. (2016)]. Beside this emerging chaotic behavior in our simulations other phenomena occurred in the parameter space exploration: a state of hypersynchronization between regions (Figures 5H,J) and a state of a “zero-line” with no oscillations that clearly does not reflect a physiological brain state (Figures 5I,J).

In order to locate the individual simulations in the spectrum of possible dynamics, meaning in the range of possible Abeta load, we examined extreme values of Abeta distribution. The virtual brains with a mean Abeta load of zero (Supplementary Figure 3A) and with the maximum Abeta load at all regions (Supplementary Figure 3B), we see as expected for the Abeta-free system a behavior similar to the low-Abeta-containing HC participants. This is not surprising, because when the HC subjects do not have a high Abeta signal, the dynamics will converge to those with zero Abeta, which is in fact then only determined by the underlying standard SC and therefore remains the same for all participants. However, the homogeneous application of maximum Abeta burden does not lead to an AD-like pattern but shows a zero-line at the whole spectrum.

To give a mathematical explanation of those phenomena, we related each participants Abeta-burden to the corresponding inhibitory time constant τ_i and used former analyses of the uncoupled local Jansen-Rit model (Spiegler et al., 2010) to estimate the bifurcation diagrams for the coupled system in this study (Figure 6). Shown diagrams allow to predict and explain the occurrence of alpha and theta rhythms or zero-lines depending on the underlying Abeta burdens. The variation of τ_i by local Abeta burden fundamentally influences the systems bifurcations by shifting the bifurcation point along the range of external input to the pyramidal cells. As a consequence, different values of Abeta lead to a variable occurrence of two limit cycles and a stable focus. Therefore, for a single region with constant external input on pyramidal cells, depending on Abeta the region might be in an alpha limit cycle, in a theta limit cycle, in a bistable condition where both cycles are possible or in a stable focus.

FIGURE 6

Figure 6. Exemplary bifurcation diagrams of the Jansen-Rit model for three different inhibitory time constants linked to three different local Abeta burdens. The modulation of the inhibitory time constant τ_i by Abeta induces shifts in the corresponding bifurcation diagrams. All bifurcation diagrams (A,E,I) show the postsynaptic potential υ₃₀ of pyramidal cells (y) depending on the pyramidal input (x) for uncoupled simulations modified from Spiegler et al. (2010). Contrary to the implementation we used for our present TVB modeling approach, here the bifurcation diagrams explore the behavior in an uncoupled system, and in accordance with (Spiegler et al., 2010) the IPSP amplitude coefficient Hi changes inversely to τ_i to keep the product of synaptic gains and dendritic time constants constant. The default input m_3T,0 on pyramidal cells starts at a firing rate of 108.5/s. Because of the potential-to-firing-rate transfer function (Equation 11), global scaling factor G is affecting both the input currents and the firing rates. For higher values of G, the input on pyramidal cells is expected to increase. First Columns, panels (A–D): Bifurcation diagram with the default time constant of 14 ms. This appears in the simulation if the Abeta SUVR is below the clinical cut-off 1.4, because then the time constant is unaffected according the transfer function in Equation 14. An Abeta burden below 1.4 SUVR corresponds to THAL phase 0 (Murray et al., 2015), i.e., that we expect no Abeta pathology in such a brain region. In this situation, there is only one limit cycle existing, which produces a frequency in alpha range (A). After increasing the input on the pyramidal cells, the alpha cycle collapses and transforms to a stable focus, where no oscillations appear in the absence of noise. This is the “zero-line” in our results. (B,D): HC participant 22 shows monomorphic alpha for lower G (green and blue line) and zero-line for higher G (red line). The distribution of regions with this dynamical regime is shown in panel (C): almost all regions of participant 22 are in this “alpha regime” with an inhibitory time constant between 14 and 20 ms [red columns in panel (C)]. This homogeneity explains the low variance of rhythms shown in the lower G ranges of (B), because all regions are in the same limit cycle and in the absence of artificial noise there is no possibility for an amplitude modulating factor. The “alpha regime” appears for all regions with an Abeta PET SUVR below 1.95. This corresponds to neuropathologic THAL phases 0, 1, 2, and 3, i.e., the regions will have no severe Abeta pathology (Murray et al., 2015). Second column, panels (E–H): Bifurcation diagram with a time constant of 22 ms, which corresponds to an intermediate Abeta load and a bistable dynamical regime which occurs for time constants between 20 and 28 ms. This corresponds to Abeta PET SUVRs between 1.95 and 2.15 (THAL phase 4 or 5), i.e., moderate-to-severe Abeta pathology (Murray et al., 2015). (E) Starting at the blue line (initial condition in alpha cycle), with an increased input on the pyramidal cells (e.g., by the network) it gets possible to reach the second limit cycle, which produces a theta rhythm and coexists with the alpha cycle while the pyramidal input is in a specific range (120/s−170/s). When the input is increased too much (e.g., by many connections of the network or by increased coupling factor G), the theta cycle disappears and the system jumps back to the alpha cycle and later on to the stable focus, which shows no oscillations in the absence of noise. This can explain some of the spectral behaviors we observed typically in the AD group (F,H): It starts with chaotic rhythms in alpha (blue line) and theta (red line) and in the shown AD participant 1 then gets synchronized to either alpha or theta. With higher couplings, the frequency gets more probably synchronized to alpha (green line), because higher G indicates a higher pyramidal input and therefore a higher attraction of the alpha cycle. (G) Remarkably for the shown participant is the fact that the bistable behavior is caused by a very small amount of regions in bistable regime, which propagate the theta rhythm to most other regions in the area 200 < G < 300. Third column, panels (I–L): Bifurcation diagram with a time constant of 50 ms, which correlates to a 95th percentile Abeta load and above. Those high Abeta burdens lead to a theta dynamical regime, which occurs for time constants between 28 and 50 ms. In clinicopathology, this corresponds to Abeta PET SUVRs above 2.15—about 90% of those patients are can be classified as THAL phase 5, i.e., they have severe Abeta pathology (Murray et al., 2015). In comparison to panel (E), the alpha limit cycle disappeared in panel (I). Therefore, we expect only theta rhythms or an activity at the stable focus. The theta cycle now begins shortly above the initial condition of pyramidal input without the alpha cycle in between. For an initial input of 108.5/s the system is in a stable focus. This may explain why in the simulation with maximum Abeta load at all regions (so each with a time constant of 50 ms) we see a zero-line without alpha at lower G values (Supplementary Figure 1). (J,L) A state of theta-only rhythm appeared in few AD participants at higher Gs (blue line). In the spectral behavior of AD participant 7, we can moreover observe a strong bistable pattern with chaotic frequency distributions for G < 300. This is likely caused by the high amount of bistable regions (K), while the synchronization to theta in higher G is an effect of the high proportion of regions in theta regime.

Simulated EEG Slowing in AD Is Caused by Heterogeneous Abeta Distribution

Figure 7 displays how the mean dominant rhythms differ between the groups. In the range below G = 100 we find a slowing in the AD group. Since in the range of lower G all three groups exhibit realistic frequency spectra and no zero-lines we consider this range of G as “physiological.” Significant differences appear between AD and non-AD for ranges of high and low G and also for high alpha and low theta rhythms (Figure 7). The heterogeneous distribution of Abeta (in contrast to an averaged homogeneous distribution) plays a crucial role in the development of this AD-specific slowing. This is indicated by simulations with the mean averaged Abeta of each participant mapped on all regions. The simulations revealed a regionally more homogenous behavior in all groups (Supplementary Figure 4). Moreover, with homogeneous distribution of Abeta the slowing in AD participants does not appear: we don't see a significant change in the theta band (Figure 7B). This is a strong indicator for the importance of the individual Abeta distribution and a proof for the necessity of heterogeneous excitotoxic effects for the creation of neural slowing. However, the absence of slowing in the simulations with homogeneous Abeta distribution does not proof the importance of a specific spatial Abeta pattern. In contrast, it only shows that there must be few regions with very high Abeta pathology to slow down the system (see Figures 6F–G).

FIGURE 7

Figure 7. AD-specific slowing in EEG and LFP and influence of the heterogeneous pattern of Abeta distribution to the spectral behavior. (A,B) The panels show the “spectrograms,” more precise the amount of regions with a dominating frequency averaged for all G values and the subjects of each group. Below, black bars are indicating significant differences for all 90 examined frequencies by a Kruskal-Wallis test (compared were the means of the amount of regions in each group having this particular frequency). In (A), for the empirical Abeta distribution pattern, the red dotted line (AD) diverges from the non-AD participants with a strong presence of dominating theta (peak at 4 Hz) and the absence of zero-line rhythm (except of very few regions, see arrow). Significant differences only appear between AD and each HC and MCI, namely for high alpha/low beta and for theta/delta (black bars). At f = 1.2 Hz (red bar), the significance level is also achieved when using a strict Bonferroni correction (p < 0.05/90). In contrast, (B) shows the same plot if the spatial distribution was “blurred”: There is no visual difference between the behavior of the three groups, and also no theta rhythm is existing in the simulations. All groups have a dominating zero-line behavior averaged across the full G range (see arrow). However, there are some frequencies that significantly differ between AD and each HC and MCI in alpha / beta range, which could be also visually related to small peaks at the plots beside. In theta and delta, where we would expect to see the slowing, there is no significant difference at all. Due to readability, for (A,B) the y-axis was limited to the amount of 100 regions. In (A), the zero-line peak of HC and MCI ends at 211, in panel (B) all zero-line peaks end at 323. The different spectra lead to different G-dependent mean frequencies for the groups, which significantly differ in areas of high and low G: (C,D)—comparison of EEG and LFP between groups. Mean dominant rhythms across all simulated EEG channels (C) and region-wise simulated LFPs (D) for all analyzed global coupling values. The frequencies of AD patients are significantly different in EEG as well as in the regional neuronal population signal. Filled shapes and thin lines represent the quantiles at 0.95 and 0.05 for each group. (C) For EEG one can see that the 95%-quantile of AD and HC as well as MCI is not overlapping in the physiological area of lower G, where AD tends to slower frequencies. In a Kruskal-Wallis test, the difference between the means of all channel frequencies per subject in the three groups is significant for AD and non-AD at 0 < G < 60 (each AD to HC and AD to MCI: p < 0.0001). They are also significantly different in the area of higher G, where AD is faster—at 450 < G < 470 (each AD to HC and AD to MCI: p < 0.0001). (D) For simulated regional neural signal the slowing effect is less prominent. The broader range of frequencies for AD is represented by the high and low limit of the 95%-quantile. This can be related to the two frequency clusters in AD at alpha and theta, which are not frequently apparent in non-AD (as in Figure 5). In a Kruskal-Wallis test, the difference between the means of all regional frequencies per subject in the three groups is only continuously significant for AD against HC at 400 < G < 450 (AD compared to HC p < 0.0001). For the other comparisons, only isolated G values deliver significant differences in the area of low G (HC and MCI) and intermediate G (AD and MCI). Because of the big amount of tests necessary to test all global coupling values, none of the tested G values achieved Bonferroni corrected significance. However, because we assume that neither the frequencies at (A,B) nor the G values at (C,D) are independent variables (which is also the reason for the “grouped” clusters of significance at alpha and theta and G = 50 and G = 450), a Bonferroni correction is not necessary.

Intra-individual Ratio of High vs. Low Abeta Burden Across All Regions Determines Simulated EEG Frequency Spectrum—Distinct Spatial Configurations of Abeta Do Not Matter for Slowing

We next examined how LFP/EEG slowing is related to the underlying Abeta burden (Figure 8). We revealed significant linear dependencies for all groups between Abeta burden and frequency. We found a strong inverse dependency for AD (R² = 0.625), i.e., an Abeta-dependent EEG slowing. In contrast, for non-AD participants the relation was revers, i.e., higher values of Abeta caused EEG acceleration.

FIGURE 8

Figure 8. Abeta-dependent slowing of LFPs is specific for AD participants. Meanwhile there is a significant linear dependency between Abeta and LFP frequency for all groups, only for AD a higher burden of Abeta leads to a decrease of frequency. HC and MCI show inverse correlations. Plotted are density plots showing the dependency between the local Abeta loads and LFPs. (A) HC group, (B) MCI group, and (C) AD group. The matrices are containing the resulting regional peak frequencies for all examined coupling values G for all participants. Linear regressions (black lines) revealed highly significant regression coefficients (p < 0.0001). A strong linear dependency between mean Abeta and LFP, that explains the greater part of the variance, is only apparent in the AD group (C). 37.5% of the variance yet cannot be explained by this linear dependency. Moreover, only for AD the dependency leads to slower frequencies for higher Abeta SUVRs, meanwhile HC and MCI have slightly faster frequencies for higher Abeta SUVRs. Visually one can see at least four contributing patterns in the AD group (C): (1) the linear decrement of frequency for higher Abeta, shown by the regression line, (2) the two frequency clusters (orange spots) at alpha and theta, (3) some regions with the zero-line behavior, particular those with low Abeta (thin line at the left, with SUVR of about 1.5), and (4) a broad variability of frequencies for regions of the same Abeta SUVR (horizontal distribution). These phenomena cannot be explained completely by a linear dependency and moreover not by a linear system at all. The criticality that divides the dynamics into three different frequency modes (zero, alpha and theta) is a phenomenon of the Jansen-Rit model as a non-linear system (Figure 6 and Supplementary Figure 6) and the broad frequency distribution is (probably) a network effect.

To test if specific regions are more important for the observed phenomena, we had to overcome the bias that only specific regions were strongly affected by Abeta. I.e., for the empirical Abeta distribution we cannot say e.g., for a region with high Abeta if it shows EEG/LFP slowing only because of its high Abeta value or because of its specific spatial and graph theoretical position in the network. Therefore, we next performed simulation with 10 random spatial distributions of the individual Abeta PET SUVRs for the 10 AD participants. In these simulations, the neural slowing appeared similarly to the empirical spatial distributions of Abeta (Supplementary Figure 5), which indicates a minor role of the distinct spatial patterns of Abeta. Instead, the ratio of regions corresponding to the three different dynamical regimes (alpha, theta, and bistable) determined the simulated frequency spectrum (Supplementary Figure 6). For an optimal value of G with 100 < G < 150, the ratio of regions with an Abeta value in theta regime best corresponded to the ratio of regions with theta frequency in LFP. Moreover, the number of regions in different regimes enables to predict the individual spectral behavior across G. This analysis shows the crucial role of G for the simulation dynamics. There might exist different optima of G—dependent on what phenomenon in the simulation is of interest. But for a specific phenomenon, in this case the correspondence of underlying Abeta PET to frequencies, we can find a narrow optimum of G wherein a specific behavior occurs.

The results of random spatial distribution of Abeta PET SUVRs were also used for a parameter space exploration (Figure 9). The analysis reveals that (1) alpha rhythms are only apparent for low time constants with τ_i < 30 ms, but for the full spectrum of G, more probable for lower G values; (2) relevant amounts of bistable rhythms are only apparent for 17 ms < τ_i < 39 ms and G > 120; (3) theta rhythms are present across almost the full spectra of G and τ_i, with an equal appearance across G, but with a local minimum at τ_i ≈ 18 ms, where the system is dominated by alpha and bistable rhythms. This exploration demonstrates two major insights. First, it confirms the crucial role of τ_i for the appearance of alpha or theta rhythms as we expect it out of the (non-coupled) bifurcation diagrams of Figure 6. Network effects are present (e.g., there are theta rhythms for low values of τ_i), but play a minor role here. Second, the value of G does not significantly affect the probability of theta rhythm, except of an alpha-theta shift for low τ_i < 20 ms and higher G > 160. This is caused by the coexistence of stable focus in alpha regime and theta limit cycle in theta regime for high pyramidal input (Figures 6A,I).

FIGURE 9

Figure 9. Alpha and bistable rhythms only appear in a specific part of the parameter space between G and τ_i. This parameter space exploration was done by coupled simulations and therefore includes network effects. Frequency (by color) is presented dependent on global coupling G (x) and inhibitory time constant τ_i (y). Projections to G and τ_i are shown beside the matrix plot, here the frequencies are classified into alpha rhythm (f > 8 Hz), theta rhythm (f < 5 Hz) and bistable rhythms (5 < f < 8). No relevant proportion of zero-lines appeared in the simulations. The difference to empirical EEG classes (with slightly lower borders for theta, meaning more exactly a theta/delta rhythm) are reasonable here because of the knowledge of only two different limit cycles in the examined configuration of the Jansen-Rit model (Figure 6). This is also the reason for the classification of frequencies between 5 and 8 Hz as bistable. The exploration was non-systematically performed by using all regions of random distributed Abeta SUVR values of the 10 AD participants, with 10 iterations of randomization per participant. However, except single values of τ_i, the full spectrum of τ_i could be explored. Single empty columns are filled with neighbor columns for better readability. In principle wee the an “isle” of alpha for low coupling and low time constant, while the rest of the dynamics is dominated by theta and delta. A full frequency spectrum (also green and yellow colors) is only apparent near the borders of the alpha isle in higher coupling.

Neural Slowing Propagates to Central Parts of the Network Independently of the Spatial Abeta Distribution

In the analysis of spatial distribution in relation to the organization of the underlying SC network (Figure 10), it can be seen that unless Abeta is distributed more peripherally, the Abeta-dependent effect of neural slowing is focused to central parts of the network. Even a random distribution of Abeta SUVRs leads to this effect (Figures 10E,F), indicating that this is a network effect. Probably this phenomenon is caused because the slowing effects are not only affecting the region itself, but also its local circuitry and neighbored regions. Hubs with a high degree and many close neighbors are therefore more probable of being affected by slow rhythms propagated by other regions. To relate this to empirical facts: We know from our data (Figure 10A) that Abeta is not deposited in hubs, but more in peripheral regions of the networks. This shows, however, how the consecutive pathologic slowing effect is afterwards focused to central and important parts of the networks. A weak peripheral affection of the inhibitory system therefore disturbs the full system seriously.

FIGURE 10

Figure 10. Theta rhythms affect central parts of the network independently of the spatial distribution of Abeta. (A) Abeta PET SUVR for AD participants: the distribution is diffuse along the cortex with no strong affection of subcortical hubs. This well corresponds to the neocortical stage C of Abeta distribution (Braak and Braak, 1991, 1997; Braak et al., 2006). (B) There is no linear dependency between the Abeta SUVR and the structural degree, as the graph above already indicates. In contrast to that, (C) shows the distribution of theta rhythm, computed as the proportion of each regions simulations (201 for different values of G for 10 subjects) with dominant theta rhythm (here simplified as a frequency that is below 8 Hz and not zero, so more precise the theta-delta-band). The patterns are not consistent with those of panel (A). This indicates that not the distinct region affected by Abeta is crucial, but more its local circuitry. Moreover, one can observe that regions with a higher degree often have a high appearance of theta rhythm (D) and show a linear dependency with R² = 0.183, in contrast to the distribution of Abeta (B), which hasn't shown such a dependency. This phenomenon is stable also for the random spatial distribution of Abeta SUVRs (E): Here we see even a stronger dependency (R² = 0.29) between structural degree and theta rhythm (F). This is remarkable because (unless the spatial distribution is random) the “pathologic” theta is focused on the hubs. This indicates that there must be network effects which concentrate the appearing theta to those regions with higher degree.

Virtual Therapy With the NMDA Antagonist Memantine

The former analyses have shown that Abeta-mediated simulated hyperexcitation can lead to realistic changes of simulated brain imaging signals in AD such as EEG slowing (Figures 5, 6). We therefore wanted to know if an established way to protect the brain of the hyperexcitation, which is the NMDA antagonist memantine, can lead to functional reversibility.

The idea in our model is now that in theory memantine acts anti-excitotoxic via its NMDA antagonism and should therefore be able to weaken the hyperexcitation we introduced to the system by Abeta (Figure 11). As mentioned above, the local coupling parameter c₃₁ represents the main part of the glutamatergic transmission and can therefore also be seen as a surrogate of NMDAergic transmission (Figure 11A). We homogeneously increased the default value of c₃₁ stepwise to observe the effects on the system. The analysis of the Jansen-Rit equilibria supports the concept of lower excitation introduced by lower c₃₁ (Figures 11B,C). The equilibrium manifold is the manifold of fixed points (the equilibrium) that is projected onto the PSP at the pyramidal cells as a function of two parameters, that is, the local excitatory-to-pyramidal coupling coefficient c₃₁ and the input on pyramidal cells m_{3T, 0}. The manifold is the object onto which the system is moving or collapsing dependent on the parameters—in a way the equilibrium that underlies the dynamics of the system. The virtual memantine leads to a partial reversibility of the altered dominant frequencies in AD compared to HC/MCI (Figure 11D). Virtual memantine increases the mean dominant EEG frequency. These simulated functional effects do not imply reversibility of neurodegeneration, but they illustrate how pharmacological intervention can theoretically counteract those processes. This observation provides first a potential mechanistic explanation of the pharmacodynamics of memantine. Second, it shows that TVB in general and the Abeta-hyperexcitation model of this study in particular are able to test the efficacy of treatment strategies such as drugs and have therefore the potential to be used for the discovery of new treatment options. Finally, it supports the concept of this study, where the impaired inhibitory function is modeled by an increased synaptic delay and it indirectly indicates that higher Abeta (by increasing τ_i) has led to a local hyperexcitation. It is to mention, that in an uncoupled network both the decrease of c₃₁ (memantine) and the increase of τ_i (by Abeta) would have the same effect (Figures 11B,C, 13C,D). In a coupled simulation, the effects are in contrast antagonistic. One reason for this seems to be, that the effect of virtual memantine is focused to central parts of the network (Figures 11E,F)—the same parts, where slowing (Figure 10) and Abeta-derived hyperexcitation (Figure 12) are occurring. The homogeneously applied memantine evolves its action, guided by the topology of the SC network, along the same path as the hyperexcitation is distributed. The effects of altered delay of GABA transmission can be reversed by adjusting NMDA transmission at another subset of the local population model. This illustrates that theoretically an alteration of the inhibitory transmission dynamics may lead to disinhibition causing hyperexcitation in downstream populations, which is reversible by reduction of excitatory input.

FIGURE 11

Figure 11. Modeling NMDA antagonism by virtual memantine. We modified the local dynamics for the AD group by homogeneously decreasing the coefficient c₃₁, which represents the coupling from excitatory population to the pyramidal cell and therefore is a potential surrogate for NMDA receptor activity (A). The coefficient c₃₁ was decreased by 25% to model the effect of memantine and was applied homogenously to all regions for the 10 AD participants. (B) Equilibrium manifold of the Jansen-Rit model, wherein υ₃₀ is a function of the model input m_{3T, 0} and c₃₁ for the median Abeta load of AD participants β = 2.1447. Note the decrease of PSP at pyramidal cells υ₃₀ with increasing c₃₁ for a constant input level—this can also be seen in the top view of the same three-dimensional plot in panel (C), where our view is parallel to the z axis. The slope of the manifold with the input decreases with increasing c₃₁. Blue areas indicate the lower branch of the equilibrium manifold and red areas the upper branch (white areas are unstable). This demonstrates, as we suggested, that when maintaining the same input level, a lower c₃₁ leads to a lower PSP at the pyramidal cells. (D) Mean EEG frequency for the three groups HC (blue), MCI (green), AD (red, with shadowed area for the range between 5th and 95th percentile) and AD with memantine (red dotted line). The virtual application of memantine shifts the AD group to the level of HC and MCI (arrow) and out of the variance of AD without memantine. (E) Change in the relative PSP (local PSP relative to mean of all regions) due to virtual memantine in AD patients. The decrement of activity seems to be focused to the hubs and to central parts of the network. Interestingly, this is the same topological area wherein relative Abeta-derived hyperactivity takes place (Figure 12). (F) Shows a significant and strong linear dependency (R² = 0.783) between the effect of memantine and the relative hyperactivity of regions. The homogenously applied virtual memantine therefore acts selectively in those regions, where hyperexcitation is already there. These regions in central network parts are also those, where slowing effects appears (Figure 10). The implementation of virtual memantine provides therefore a link between its supposed anti-excitotoxic effect (because it acts selectively in relatively hyperactive regions) and the reversed functional phenomenon of slowing (because slowing, hyperactivity and the memantine effect are all focused to the same parts of the network).

FIGURE 12

Figure 12. Local hyperexcitation is introduced by Abeta and spatially linked to LFP slowing. (A–D) Shows the relative firing rate of pyramidal cells, i.e., the difference between the local firing rate and the global firing rate (averaged across all regions for one simulation / one value of G). Because the maximum firing rate is 5/s, the relative firing rate ranges from −5/s to 5/s. The firing rate can be calculated by using the potential-to-firing-rate transfer function from Equation 11. Color indicates the natural logarithm of the number of regions at each point of the histogram. (A) In the absence of Abeta, where all regions have an inhibitory time constant τ_i = 14 ms, the firing rate shows low heterogeneity. There are neither hypoactive, nor hyperactive regions—the whole systems activity is near the “baseline of the brain” (mean firing rate of all regions). (B) For HC participants, the Abeta burden of some regions already introduces heterogeneity. Although, most regions are still near the baseline, as also in (C) for the MCI participants. In panel (D), we can see that the AD participants have enough Abeta to introduce a strong local heterogeneity of firing rate. There are now about 4 clusters, which are stable across G: one near the baseline [but in contrast to (A–C) not directly at the baseline, but below], another hypoactive and two hyperactive clusters of regions. The spatial distribution of the relative activity is demonstrated in (E): Here we can observe that the hyperexcitation of pyramidal cells is focused to the hubs and to central parts of the network, which is also shown in (G). This pattern is similar to what we could observe for the theta rhythm, which is also focused to central network parts (Figure 10) and to the effect of memantine (Figure 11). There is moreover a significant linear dependency between the relative PSP and the probability of theta rhythms (F) as well as a significant and strong linear dependency (R² = 0.594) between the relative PSP and the natural logarithm of the structural degree (G).

FIGURE 13

Figure 13. Evaluation of the used cause-and-effect model. (A) Inhibitory convolution kernel h_i(t) as a function of Abeta. The kernel is flattened with increasing Abeta and the area under the curve increases as follows: AUC = H_i • τ_i(β). Therefore, longer τ_i by higher Abeta leads to a slowed down filter action (because of the delayed maximum of the IPSP). As a side effect, because of the higher AUC in the inhibitory transmission, the overall inhibition in the system increases—this can also be seen in (C) and (D): Equilibrium manifold (the manifold of fixed points onto which the system is collapsing) of the Jansen-Rit model in three-dimensional perspective (C) and in top view in parallel to the z axis (D), wherein υ₃₀ is a function of the model input m_{3T, 0} and β. The reaction of the system to input is changed by Abeta. Note the decrease of PSP at pyramidal cells υ₃₀ with increasing β for a constant input level—this can also be seen in (C). The slope of the manifold with the input decreases with increasing β. Blue areas indicate the lower branch of the equilibrium manifold and red areas the upper branch (white areas are unstable). This demonstrates, as we suggested because of (A), that when maintaining the same input level, a higher Abeta burden leads to a globally lower PSP baseline at the pyramidal cells. Because of the larger AUC at inhibitory transmission, we need more input to achieve the same global activity. As we would not expect the brain to increase or decrease its global activity level and are more interested in the local activity, this can be seen as a mathematical artifact as well as a limitation of the model. To evaluate local hyperexcitation, we therefore used the relative PSP to the baseline in Figure 12 instead of the absolute PSP. (B) The bifurcation diagram of local bifurcations shows the PSP of pyramidal cells as an explicit function of Abeta. We can see that the richest dynamic repertoire appears in the range around β = 2. Here appears multistability, which can be seen because of the folded Hopf curve. The local bifurcations support our hypothesis of diversified regional activity introduced by Abeta: above β = 2, two stable foci coexist, while one of them is positioned at a low PSP and the other one at a high PSP. This enables us to expect a diversified activity to nodes with lower and nodes with higher activity—as well, because the diagram shows only the local bifucrations, it does not proof it. (E) This leads together with our other results (in particular Figures 10–12) to the following mechanistic theory: longer inhibitory time constants, introduced by Abeta, lead to relatively more hypoactive and more hyperactive regions because of a disturbance of the local E/I balance in the vicinity of Abeta (Figures 12A–D). The hyperactivity is focused to the hubs and central network parts (Figures 12E–G). Because we know that this pattern does not evolve from the spatial Abeta distribution (Figures 10A,B), we can say that the hyperactivity is centralized (to central network parts) and the hypoactivity is vice versa peripherialized (to the periphery of the graph). This phenomenon goes along with the distribution of slow rhythms in our results (Figure 10). Therefore, our model linked local disturbed E/I balance with a slowed filter action of the inhibitory signal to hyperactivity along with slowing in central parts of the brain network. Parts of the figure are modified from Deco et al. (2017).

Discussion

Local Abeta-mediated disinhibition and hyperexcitation are considered candidate mechanisms of AD pathogenesis. In TVB simulations, the molecular candidate mechanism has led to macro-scale slowing in EEG and neural signal with a particular shift form alpha to theta previously observed in AD patients (Loewenstein et al., 1989; Rice et al., 1990; Dierks et al., 2000; Babiloni et al., 2016; Malek et al., 2017; Nakamura et al., 2018). These observations cannot be directly inferred by the hyperexcitation implemented in our model. Because we standardized all other factors and used a common SC for all simulations, this approach enables to examine the effects of locally altered E/I balance on an individual whole-brain level but without any other confounding factors.

We showed that the slowing in simulated EEG and LFP is specific for the AD group (Figures 7, 8). This offers an explanation, how the shift from alpha to theta, that is observable in EEG of AD patients (Loewenstein et al., 1989; Rice et al., 1990; Dierks et al., 2000; Babiloni et al., 2016; Malek et al., 2017; Nakamura et al., 2018), could be explained on a synaptic level—namely by an impaired inhibition. This computational modeling result supports the findings of specific toxicity of Abeta to inhibitory neurons (Ripoli et al., 2014; Ulrich, 2015; Ren et al., 2018).

We demonstrate the computational principles underlying this Abeta dependent slowing of EEG/LFP (Figure 6, Supplementary Figure 6). Dependent on the Abeta burden alpha, theta or bistable regime develop caused by an alteration of the inhibitory time constant that leads to changes of the systems bifurcation behavior (Figure 6, Supplementary Figures 1, 2, 6).

The simulated LFP/EEG slowing in AD patients crucially depends on the spatially heterogenous Abeta distribution as measured by PET—the slowing disappears when using a homogenously distributed mean Abeta burden instead for simulation (Figure 7). To exhibit the slowing effect few regions with high Abeta burden are required—while the specific location of these regions seems not to be relevant (Supplementary Figure 5). The crucial factor for AD-specific slowing behavior in our simulations is the presence of very few regions that are strongly affected by Abeta (Figure 6F and middle column of Supplementary Figure 6).

Independently of the location of high Abeta burdens in the simulated brain, slowing emerges at the core, i.e., hubs of the structural connectome (Figure 10). This indicates that the central parts of the system are impacted functionally by the Abeta burden. Moreover, it shows that while Abeta is often distributed in peripheral parts of the structural connectome, its functional consequences affect the important hubs. This could provide a possible explanation why a peripheral distribution of Abeta leads to severe disturbances of cognitive function.

Abeta leads by the disturbance of E/I balance to more local hyperexcitation (Figure 12). Because the range of activity is broader, we have more hypoactive and also more hyperactive regions (Figure 13). The local hyperexcitation is strongly corrleated with local LFP slowing (Figure 10) and also focused to the hubs of the network (Figure 12).

We also showed that the drug memantine that is known for improving brain function in severe AD can be modeled by a decreased transmission between the excitatory interneurons and the pyramidal cells and is able to achieve a “normalized” brain function in silico, too (Figure 11). Its effect is evolves selectively in hyperactive regions and in those parts of the network, where slow rhythms appear. This moreover demonstrates the potential of TVB to test and develop new treatment strategies.

One major limitation of this study is the lack of direct validation of the simulated electrophysiological phenomena. Neither EEG nor LFP data was available in the ADNI-3 cohort. Although EEG slowing in AD is an established concept (Loewenstein et al., 1989; Rice et al., 1990; Dierks et al., 2000; Babiloni et al., 2016; Malek et al., 2017; Nakamura et al., 2018), future studies will have to validate the simulated data directly with individual EEG.

The second important limitation is the implementation of disturbed E/I balance by the inhibitory time constant. Although the longer time constants lead to slowed filter action (Figure 13A) and local hyperexcitation at important network structures (Figure 12), the global activity is decreased for the same input level. To overcome this model limitation, it would be necessary to correct the input level, e.g., by increasing the default input m_{3T, 0} with higher mean Abeta burdens or by increasing a coefficient inside the IRF (Equation 2) to keep the AUC and therefore the transmitted energy at inhibitory transmission constant. This should be examined in future studies to evaluate the effect of such a correction. However, this would only be necessary if the global activity level would be a target of interest for another research question. Because of the feedback loops in a coupled brain network, each introduction of over- or dis-inhibition will lead both to hypo- and hyper-active regions. An analysis of hyperactivity needs therefore always a control activity—because hyperactivity can be meant spatially, temporally, or dependent on other factors. In our model, we could introduce spatially distributed hyperactivity (Figure 12) that was linked to local slowing (Figure 10), network topology (Figures 12E–G) and could be antagonized selectively by virtual memantine (Figure 11).

The differential equations that describeThe differential equations that describeOf course the pathophysiological model presented in this study can only cover a small subspace of possible AD mechanisms. Even Abeta itself is probably only one player in the multifactorial pathogenesis (Selkoe and Hardy, 2016; Gauthier et al., 2018), and synaptic disinhibition is only one of its effects (Furukawa et al., 1996; Good et al., 1996; Hardy and Selkoe, 2002; Chen, 2005; Busche et al., 2012; Prasansuklab and Tencomnao, 2013; Bloom, 2014; Sadigh-Eteghad et al., 2015; Song et al., 2015; Ulrich, 2015; Celebi et al., 2016; Gauthier et al., 2018; Ren et al., 2018). We are fully aware that other major hallmarks as Tau (Bloom, 2014; Guo et al., 2017; Tapia-Rojas et al., 2019) and cardiovascular risk factors (Love and Miners, 2016; Storck and Pietrzik Claus, 2018; Bannai et al., 2019) cannot be excluded in the discussion of AD etiology, as well as alternative concepts such as microglia and neuroinflammation (Heneka et al., 2015a,b; Wang and Colonna, 2019; Zhou et al., 2019), polygenetic risk factors (Mahley, 2016; Hudry et al., 2019; Jansen et al., 2019; Takatori et al., 2019), environmental factors as neurotoxic or infectious agents (Alonso et al., 2018; McLachlan et al., 2019), and concomitant proteinopathies (Robinson et al., 2018a,b).

Another limitation is the small sample size of 33 participants. Future studies will have to consider much more participants, which will then help to formulate even more general conclusions. However, because of emergent effects in the brain simulation, differences between the groups were often very relevant and significant. An overview of all used statistical test in this study can be found in Supplementary Table 8.

However, we present a first proof of concept for linking molecular changes as detected by PET to large-scale brain modeling using the simulation framework TVB. This study therefore can work as a blueprint for future approaches in computational brain modeling bridging scales of neural function. For the research on AD pathogenesis, this study provides a possible mechanistic explanation that links Abeta-related synaptic disinhibition at the micro-scale to AD-specific EEG slowing. In general, our study can be seen as proof of concept that TVB enables research on disease mechanisms at a multi-scale level and has potential to lead to improved diagnostics and to the discovery of new treatments.

Data Availability

The raw data for this study is available in ADNI. The codes used in this study are available on request to the corresponding author.

Ethics Statement

This study has been approved from the Ethics Board of the Charité - Universitätsmedizin Berlin under the approval number EA2/100/19.

Author Contributions

All authors have made substantial intellectual contributions to this work and approved it for publication. LS and PR had the idea to this study. LS, PT, ASp, and PR developed the concept and study design. LS wrote the manuscript, conducted the analysis and interpretation of results, and developed the figures. PT performed the MRI and PET image processing and supercomputer simulations. PT, ASp, M-AD-C, ASo, VJ, AM, and PR contributed to the interpretation of the results, figure development, and writing of the manuscript.

Funding

Data collection and sharing for this project was funded by the Alzheimer's Disease Neuroimaging Initiative (ADNI) (National Institutes of Health Grant U01 AG024904) and DOD ADNI (Department of Defense award number W81XWH-12-2-0012). ADNI is funded by the National Institute on Aging, the National Institute of Biomedical Imaging and Bioengineering, and through generous contributions from the following: AbbVie, Alzheimer's Association; Alzheimer's Drug Discovery Foundation; Araclon Biotech; BioClinica, Inc.; Biogen; Bristol-Myers Squibb Company; CereSpir, Inc.; Cogstate; Eisai Inc.; Elan Pharmaceuticals, Inc.; Eli Lilly and Company; EuroImmun; F. Hoffmann-La Roche Ltd and its affiliated company Genentech, Inc.; Fujirebio; GE Healthcare; IXICO Ltd.; Janssen Alzheimer Immunotherapy Research and Development, LLC.; Johnson and Johnson Pharmaceutical Research and Development LLC.; Lumosity; Lundbeck; Merck and Co., Inc.; Meso Scale Diagnostics, LLC.; NeuroRx Research; Neurotrack Technologies; Novartis Pharmaceuticals Corporation; Pfizer Inc.; Piramal Imaging; Servier; Takeda Pharmaceutical Company; and Transition Therapeutics. The Canadian Institutes of Health Research is providing funds to support ADNI clinical sites in Canada. Private sector contributions are facilitated by the Foundation for the National Institutes of Health (www.fnih.org). The grantee organization is the Northern California Institute for Research and Education, and the study is coordinated by the Alzheimer's Therapeutic Research Institute at the University of Southern California. ADNI data are disseminated by the Laboratory for Neuro Imaging at the University of Southern California.

PR acknowledges the following funding sources: H2020 Research and Innovation Action grants 826421 and 650003, 720270 and 785907, and ERC 683049; German Research Foundation CRC 1315 and 936, and RI 2073/6-1; Berlin Institute of Health and Foundation Charité, Johanna Quandt Excellence Initiative.

We acknowledge support from the German Research Foundation (DFG) and the Open Access Publication Fund of Charité—Universitätsmedizin Berlin.

Further we acknowledge Lea Doppelbauer and Jan Roediger for their helpful discussions.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fncom.2019.00054/full#supplementary-material

References

Alonso, R., Pisa, D., Fernandez-Fernandez, A. M., and Carrasco, L. (2018). Infection of fungi and bacteria in brain tissue from elderly persons and patients with Alzheimer's disease. Front. Aging Neurosci. 10:159. doi: 10.3389/fnagi.2018.00159

PubMed Abstract | CrossRef Full Text | Google Scholar

Alzheimer's Association (2018). 2018 Alzheimer's disease facts and figures. Alzheimer's Dement. 14, 367–429. doi: 10.1016/j.jalz.2018.02.001

CrossRef Full Text | Google Scholar

Babiloni, C., Del Percio, C., Caroli, A., Salvatore, E., Nicolai, E., Marzano, N., et al. (2016). Cortical sources of resting state EEG rhythms are related to brain hypometabolism in subjects with Alzheimer's disease: an EEG-PET study. Neurobiol. Aging 48, 122–134. doi: 10.1016/j.neurobiolaging.2016.08.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Bannai, T., Mano, T., Chen, X., Ohtomo, G., Ohtomo, R., Tsuchida, T., et al. (2019). Chronic cerebral hypoperfusion shifts the equilibrium of amyloid beta oligomers to aggregation-prone species with higher molecular weight. Sci. Rep. 9:2827. doi: 10.1038/s41598-019-39494-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Biogen (2019). Biogen and EISAI to discontinue phase 3 engage and emerge trials of aducanumab in Alzheimer's disease 2019 [04/01/2019]. Available online at: http://investors.biogen.com/news-releases/news-release-details/biogen-and-eisai-discontinue-phase-3-engage-and-emerge-trials

Blennow, K., de Leon, M. J., and Zetterberg, H. (2006). Alzheimer's disease. Lancet 368, 387–403. doi: 10.1016/S0140-6736(06)69113-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Bloom, G. S. (2014). Amyloid-β and tau: The trigger and bullet in Alzheimer disease pathogenesis. JAMA Neurol. 71, 505–508. doi: 10.1001/jamaneurol.2013.5847

PubMed Abstract | CrossRef Full Text | Google Scholar

Bojak, I., Oostendorp, T. F., Reid, A. T., and Kotter, R. (2010). Connecting mean field models of neural activity to EEG and fMRI data. Brain Topogr. 23, 139–149. doi: 10.1007/s10548-010-0140-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Braak, H., Alafuzoff, I., Arzberger, T., Kretzschmar, H., and Del Tredici, K. (2006). Staging of Alzheimer disease-associated neurofibrillary pathology using paraffin sections and immunocytochemistry. Acta Neuropathol. 112, 389–404. doi: 10.1007/s00401-006-0127-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Braak, H., and Braak, E. (1991). Neuropathological stageing of Alzheimer-related changes. Acta Neuropathol. 82, 239–259. doi: 10.1007/BF00308809

PubMed Abstract | CrossRef Full Text | Google Scholar

Braak, H., and Braak, E. (1997). Frequency of stages of Alzheimer-related lesions in different age categories. Neurobiol. Aging 18, 351–357. doi: 10.1016/S0197-4580(97)00056-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Breakspear, M. (2017). Dynamic models of large-scale brain activity. Nat. Neurosci. 20, 340–352. doi: 10.1038/nn.4497

PubMed Abstract | CrossRef Full Text | Google Scholar

Busche, M. A., Chen, X., Henning, H. A., Reichwald, J., Staufenbiel, M., Sakmann, B., et al. (2012). Critical role of soluble amyloid-beta for early hippocampal hyperactivity in a mouse model of Alzheimer's disease. Proc. Natl. Acad. Sci. U.S.A. 109, 8740–8745. doi: 10.1073/pnas.1206171109

PubMed Abstract | CrossRef Full Text | Google Scholar

Busche, M. A., Eichhoff, G., Adelsberger, H., Abramowski, D., Wiederhold, K. H., Haass, C., et al. (2008). Clusters of hyperactive neurons near amyloid plaques in a mouse model of Alzheimer's disease. Science 321, 1686–1689. doi: 10.1126/science.1162844

PubMed Abstract | CrossRef Full Text | Google Scholar

Cabral, J., Kringelbach, M. L., and Deco, G. (2017). Functional connectivity dynamically evolves on multiple time-scales over a static structural connectome: models and mechanisms. Neuroimage 160, 84–96. doi: 10.1016/j.neuroimage.2017.03.045

PubMed Abstract | CrossRef Full Text | Google Scholar

Celebi, O., Uzdogan, A., Oguz, K. K., Has, A. C., Dolgun, A., Cakmakli, G. Y., et al. (2016). Default mode network connectivity is linked to cognitive functioning and CSF Abeta1-42 levels in Alzheimer's disease. Arch. Gerontol. Geriatr. 62, 125–132. doi: 10.1016/j.archger.2015.09.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, C. (2005). beta-Amyloid increases dendritic Ca2+ influx by inhibiting the A-type K+ current in hippocampal CA1 pyramidal neurons. Biochem. Biophys. Res. Commun. 338, 1913–1919. doi: 10.1016/j.bbrc.2005.10.169

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, R., Chan, P. T., Chu, H., Lin, Y. C., Chang, P. C., Chen, C. Y., et al. (2017). Treatment effects between monotherapy of donepezil versus combination with memantine for Alzheimer disease: a meta-analysis. PLoS ONE 12:e0183586. doi: 10.1371/journal.pone.0183586

PubMed Abstract | CrossRef Full Text | Google Scholar

Chiao, P., Bedell, B. J., Avants, B., Zijdenbos, A. P., Grand'Maison, M., O'Neill, P., et al. (2019). Impact of reference and target region selection on amyloid PET SUV ratios in the phase 1b PRIME study of Aducanumab. J. Nucl. Med. 60, 100–106. doi: 10.2967/jnumed.118.209130

PubMed Abstract | CrossRef Full Text | Google Scholar

Clark, C. M., Schneider, J. A., Bedell, B. J., Beach, T. G., Bilker, W. B., Mintun, M. A., et al. (2011). Use of florbetapir-pet for imaging β-amyloid pathology. JAMA. 305, 275–283. doi: 10.1001/jama.2010.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

David, O., and Friston, K. J. (2003). A neural mass model for MEG/EEG:: coupling and neuronal dynamics. NeuroImage 20, 1743–1755. doi: 10.1016/j.neuroimage.2003.07.015

PubMed Abstract | CrossRef Full Text | Google Scholar

David, O., Kiebel, S. J., Harrison, L. M., Mattout, J., Kilner, J. M., and Friston, K. J. (2006). Dynamic causal modeling of evoked responses in EEG and MEG. NeuroImage 30, 1255–1272. doi: 10.1016/j.neuroimage.2005.10.045

PubMed Abstract | CrossRef Full Text | Google Scholar

Deco, G., Jirsa, V. K., Robinson, P. A., Breakspear, M., and Friston, K. (2008). The dynamic brain: from spiking neurons to neural masses and cortical fields. PLoS Comput. Biol. 4:e1000092. doi: 10.1371/journal.pcbi.1000092

PubMed Abstract | CrossRef Full Text | Google Scholar

Deco, G., Kringelbach, M. L., Jirsa, V. K., and Ritter, P. (2017). The dynamics of resting fluctuations in the brain: metastability and its dynamical cortical core. Sci. Rep. 7:3095. doi: 10.1038/s41598-017-03073-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Dierks, T., Jelic, V., Pascual-Marqui, R. D., Wahlund, L., Julin, P., Linden, D. E., et al. (2000). Spatial pattern of cerebral glucose metabolism (PET) correlates with localization of intracerebral EEG-generators in Alzheimer's disease. Clin. Neurophysiol. 111, 1817–1824. doi: 10.1016/S1388-2457(00)00427-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Doggrell, S. A. (2018). Grasping at straws: the failure of solanezumab to modify mild Alzheimer's disease. Expert Opin. Biol. Ther. 18, 1189–1192. doi: 10.1080/14712598.2018.1543397

PubMed Abstract | CrossRef Full Text | Google Scholar

Falcon, M. I., Riley, J. D., Jirsa, V., McIntosh, A. R., Elinor Chen, E., and Solodkin, A. (2016). Functional mechanisms of recovery after chronic stroke: modeling with the virtual brain. eNeuro 3:ENEURO.0158-15.2016. doi: 10.1523/ENEURO.0158-15.2016

PubMed Abstract | CrossRef Full Text | Google Scholar

Falcon, M. I., Riley, J. D., Jirsa, V., McIntosh, A. R., Shereen, A. D., Chen, E. E., et al. (2015). The virtual brain: modeling biological correlates of recovery after chronic stroke. Front. Neurol. 6:228. doi: 10.3389/fneur.2015.00228

PubMed Abstract | CrossRef Full Text | Google Scholar

Freeman, W. J. (ed). (1975). Mass Action in the Nervous System. San Diego, CA: Academic Press.

Google Scholar

Fruchterman, T. M. J., and Reingold, E. M. (1991). Graph drawing by force-directed placement. Softw. Pract. Exper. 21, 1129–1164. doi: 10.1002/spe.4380211102

CrossRef Full Text | Google Scholar

Furukawa, K., Barger, S. W., Blalock, E. M., and Mattson, M. P. (1996). Activation of K+ channels and suppression of neuronal activity by secreted beta-amyloid-precursor protein. Nature 379, 74–78. doi: 10.1038/379074a0

PubMed Abstract | CrossRef Full Text | Google Scholar

Gauthier, S., Zhang, H., Ng, K. P., Pascoal, T. A., and Rosa-Neto, P. (2018). Impact of the biological definition of Alzheimer's disease using amyloid, tau and neurodegeneration (ATN): what about the role of vascular changes, inflammation, Lewy body pathology? Transl. Neurodegen. 7:12. doi: 10.1186/s40035-018-0117-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghosh, A., Rho, Y., McIntosh, A. R., Kotter, R., and Jirsa, V. K. (2008). Noise during rest enables the exploration of the brain's dynamic repertoire. PLoS Comput. Biol. 4:e1000196. doi: 10.1371/journal.pcbi.1000196

PubMed Abstract | CrossRef Full Text | Google Scholar

Glasser, M. F., Coalson, T. S., Robinson, E. C., Hacker, C. D., Harwell, J., Yacoub, E., et al. (2016). A multi-modal parcellation of human cerebral cortex. Nature 536, 171–178. doi: 10.1038/nature18933

PubMed Abstract | CrossRef Full Text | Google Scholar

Glasser, M. F., Sotiropoulos, S. N., Wilson, J. A., Coalson, T. S., Fischl, B., Andersson, J. L., et al. (2013). The minimal preprocessing pipelines for the human connectome project. Neuroimage 80, 105–124. doi: 10.1016/j.neuroimage.2013.04.127

PubMed Abstract | CrossRef Full Text | Google Scholar

Good, T. A., Smith, D. O., and Murphy, R. M. (1996). Beta-amyloid peptide blocks the fast-inactivating K+ current in rat hippocampal neurons. Biophys. J. 70, 296–304. doi: 10.1016/S0006-3495(96)79570-X

PubMed Abstract | CrossRef Full Text | Google Scholar

Grienberger, C., Rochefort, N. L., Adelsberger, H., Henning, H. A., Hill, D. N., Reichwald, J., et al. (2012). Staged decline of neuronal function in vivo in an animal model of Alzheimers disease. Nat. Commun. 3:774. doi: 10.1038/ncomms1783

CrossRef Full Text | Google Scholar

Guo, T., Noble, W., and Hanger, D. P. (2017). Roles of tau protein in health and disease. Acta Neuropathol. 133, 665–704. doi: 10.1007/s00401-017-1707-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Hardy, J., and Selkoe, D. J. (2002). The amyloid hypothesis of Alzheimer's disease: progress and problems on the road to therapeutics. Science 297, 353–356. doi: 10.1126/science.1072994

PubMed Abstract | CrossRef Full Text | Google Scholar

Hazra, A., Gu, F., Aulakh, A., Berridge, C., Eriksen, J. L., and Ziburkus, J. (2013). Inhibitory neuron and hippocampal circuit dysfunction in an aged mouse model of Alzheimer's disease. PLoS ONE 8:e64318. doi: 10.1371/journal.pone.0064318

PubMed Abstract | CrossRef Full Text | Google Scholar

Heneka, M. T., Carson, M. J., Khoury, J. E., Landreth, G. E., Brosseron, F., Feinstein, D. L., et al. (2015a). Neuroinflammation in Alzheimer's disease. Lancet Neurol. 14, 388–405. doi: 10.1016/S1474-4422(15)70016-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Heneka, M. T., Golenbock, D. T., and Latz, E. (2015b). Innate immunity in Alzheimer's disease. Nat. Immunol. 16:229. doi: 10.1038/ni.3102

PubMed Abstract | CrossRef Full Text | Google Scholar

Hindmarsh, J., and Rose, R. (1984). A model of neuronal bursting using three coupled first order differential equations. Proc. R. Soc. London B Biol. Sci. 221, 87–102. doi: 10.1098/rspb.1984.0024

PubMed Abstract | CrossRef Full Text | Google Scholar

Honey, C. J., Kotter, R., Breakspear, M., and Sporns, O. (2007). Network structure of cerebral cortex shapes functional connectivity on multiple time scales. Proc. Natl. Acad. Sci. U.S.A. 104, 10240–10245. doi: 10.1073/pnas.0701519104

PubMed Abstract | CrossRef Full Text | Google Scholar

Honig, L. S., Vellas, B., Woodward, M., Boada, M., Bullock, R., Borrie, M., et al. (2018). Trial of solanezumab for mild dementia due to Alzheimer's disease. N. Engl. J. Med. 378, 321–330. doi: 10.1056/NEJMoa1705971

PubMed Abstract | CrossRef Full Text | Google Scholar

Hudry, E., Klickstein, J., Cannavo, C., Jackson, R., Muzikansky, A., Gandhi, S., et al. (2019). Opposing roles of apolipoprotein E in aging and neurodegeneration. Life Sci. Alliance 2:e201900325. doi: 10.26508/lsa.201900325

PubMed Abstract | CrossRef Full Text | Google Scholar

Jack, C. R., Wiste, H. J., Weigand, S. D., Rocca, W. A., Knopman, D. S., Mielke, M. M., et al. (2014). Age-specific population frequencies of amyloidosis and neurodegeneration among cognitively normal people age 50-89 years: a cross-sectional study. Lancet Neurol. 13, 997–1005. doi: 10.1016/S1474-4422(14)70194-2

CrossRef Full Text | Google Scholar

Jack, C. R. Jr., Bennett, D. A., Blennow, K., Carrillo, M. C., Dunn, B., Haeberlein, S. B., et al. (2018). NIA-AA Research Framework: Toward a biological definition of Alzheimer's disease. Alzheimer's Dement. J. Alzheimer's Assoc. 14, 535–562. doi: 10.1016/j.jalz.2018.02.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Jack, C. R. Jr., Lowe, V. J., Weigand, S. D., Wiste, H. J., Senjem, M. L., Knopman, D. S., et al. (2009). Serial PIB and MRI in normal, mild cognitive impairment and Alzheimer's disease: implications for sequence of pathological events in Alzheimer's disease. Brain 132(Pt 5), 1355–1365. doi: 10.1093/brain/awp062

PubMed Abstract | CrossRef Full Text | Google Scholar

Jansen, B. H., and Rit, V. G. (1995). Electroencephalogram and visual evoked potential generation in a mathematical model of coupled cortical columns. Biol. Cybern. 73, 357–366. doi: 10.1007/BF00199471

PubMed Abstract | CrossRef Full Text | Google Scholar

Jansen, I. E., Savage, J. E., Watanabe, K., Bryois, J., Williams, D. M., Steinberg, S., et al. (2019). Genome-wide meta-analysis identifies new loci and functional pathways influencing Alzheimer's disease risk. Nat. Genet. 51, 404–413. doi: 10.1038/s41588-018-0311-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, M., Beckmann, C. F., Behrens, T. E., Woolrich, M. W., and Smith, S. M. (2012). FSL. Neuroimage 62, 782–790. doi: 10.1016/j.neuroimage.2011.09.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Jirsa, V. K. (2009). Neural field dynamics with local and global connectivity and time delay. Philos. Trans. A Math. Phys. Eng. Sci. 367, 1131–1143. doi: 10.1098/rsta.2008.0260

PubMed Abstract | CrossRef Full Text | Google Scholar

Jirsa, V. K., Jantzen, K. J., Fuchs, A., and Kelso, J. A. S. (2002). Spatiotemporal forward solution of the EEG and MEG using network modeling. IEEE Trans. Med. Imag. 21, 493–504. doi: 10.1109/TMI.2002.1009385

PubMed Abstract | CrossRef Full Text | Google Scholar

Jirsa, V. K., and Kelso, J. A. S. (2000). Spatiotemporal pattern formation in neural systems with heterogeneous connection topologies. Phys. Rev. E 62, 8462–8465. doi: 10.1103/PhysRevE.62.8462

PubMed Abstract | CrossRef Full Text | Google Scholar

Jirsa, V. K., Proix, T., Perdikis, D., Woodman, M. M., Wang, H., Gonzalez-Martinez, J., et al. (2017). The virtual epileptic patient: individualized whole-brain models of epilepsy spread. Neuroimage 145(Pt B), 377–388. doi: 10.1016/j.neuroimage.2016.04.049

PubMed Abstract | CrossRef Full Text | Google Scholar

Jirsa, V. K., Sporns, O., Breakspear, M., Deco, G., and McIntosh, A. R. (2010). Towards the virtual brain: network modeling of the intact and the damaged brain. Arch. Ital. Biol. 148, 189–205. doi: 10.4449/aib.v148i3.1223

PubMed Abstract | CrossRef Full Text | Google Scholar

Jones, E., Oliphant, E., and Peterson, P. (2001). SciPy: Open Source Scientific Tools for Python. Available online at: http://www.scipy.org/ (accessed July 30, 2019).

Kirschstein, T., and Köhling, R. (2009). What is the source of the EEG? Clin. EEG Neurosci. 40, 146–149. doi: 10.1177/155005940904000305

PubMed Abstract | CrossRef Full Text | Google Scholar

Kishi, T., Matsunaga, S., Oya, K., Nomura, I., Ikuta, T., and Iwata, N. (2017). Memantine for Alzheimer's disease: an updated systematic review and meta-analysis. J. Alzheimers Dis. 60, 401–425. doi: 10.3233/JAD-170424

PubMed Abstract | CrossRef Full Text | Google Scholar

Klunk, W. E., Price, J. C., Mathis, C. A., Tsopelas, N. D., Lopresti, B. J., Ziolko, S. K., et al. (2007). Amyloid deposition begins in the striatum of presenilin-1 mutation carriers from two unrelated pedigrees. J. Neurosci. 27, 6174–6184. doi: 10.1523/JNEUROSCI.0730-07.2007

PubMed Abstract | CrossRef Full Text | Google Scholar

Kringelbach, M. L., McIntosh, A. R., Ritter, P., Jirsa, V. K., and Deco, G. (2015). The rediscovery of slowness: exploring the timing of cognition. Trends Cogn. Sci. 19, 616–628. doi: 10.1016/j.tics.2015.07.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Kunze, T., Hunold, A., Haueisen, J., Jirsa, V., and Spiegler, A. (2016). Transcranial direct current stimulation changes resting state functional connectivity: a large-scale brain network modeling study. NeuroImage 140, 174–187. doi: 10.1016/j.neuroimage.2016.02.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Limon, A., Reyes-Ruiz, J. M., and Miledi, R. (2012). Loss of functional GABA(A) receptors in the Alzheimer diseased brain. Proc. Natl. Acad. Sci. U.S.A. 109, 10071–10076. doi: 10.1073/pnas.1204606109

PubMed Abstract | CrossRef Full Text | Google Scholar

Litvak, V., Mattout, J., Kiebel, S., Phillips, C., Henson, R., Kilner, J., et al. (2011). EEG and MEG data analysis in SPM8. Comput. Intell. Neurosci. 2011:852961. doi: 10.1155/2011/852961

PubMed Abstract | CrossRef Full Text | Google Scholar

Loewenstein, D. A., Barker, W. W., Chang, J. Y., Apicella, A., Yoshii, F., Kothari, P., et al. (1989). Predominant left hemisphere metabolic dysfunction in dementia. Arch. Neurol. 46, 146–152. doi: 10.1001/archneur.1989.00520380046012

PubMed Abstract | CrossRef Full Text | Google Scholar

Logovinsky, V., Satlin, A., Lai, R., Swanson, C., Kaplow, J., Osswald, G., et al. (2016). Safety and tolerability of BAN2401 - a clinical study in Alzheimer's disease with a protofibril selective Aβ antibody. Alzheimer's Res. Ther. 8:14. doi: 10.1186/s13195-016-0181-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Lopes da Silva, F. H., Hoeks, A., Smits, H., and Zetterberg, L. H. (1974). Model of brain rhythmic activity. The alpha-rhythm of the thalamus. Kybernetik 15, 27–37. doi: 10.1007/BF00270757

PubMed Abstract | CrossRef Full Text | Google Scholar

Love, S., and Miners, J. S. (2016). Cerebrovascular disease in ageing and Alzheimer's disease. Acta Neuropathol. 131, 645–658. doi: 10.1007/s00401-015-1522-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Mahley, R. W. (2016). Apolipoprotein E: from cardiovascular disease to neurodegenerative disorders. J. Mol. Med. 94, 739–746. doi: 10.1007/s00109-016-1427-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Malek, N., Baker, M. R., Mann, C., and Greene, J. (2017). Electroencephalographic markers in dementia. Acta Neurol. Scand. 135, 388–393. doi: 10.1111/ane.12638

PubMed Abstract | CrossRef Full Text | Google Scholar

Martinez, G., Vernooij, R. W., Fuentes Padilla, P., Zamora, J., Bonfill Cosp, X., and Flicker, L. (2017). 18F PET with florbetapir for the early diagnosis of Alzheimer's disease dementia and other dementias in people with mild cognitive impairment (MCI). Cochr. Database Syst. Rev. 11:Cd012216. doi: 10.1002/14651858.CD012216.pub2

PubMed Abstract | CrossRef Full Text | Google Scholar

McKhann, G., Drachman, D., Folstein, M., Katzman, R., Price, D., and Stadlan, E. M. (1984). Clinical diagnosis of Alzheimer's disease: report of the NINCDS-ADRDA work group under the auspices of department of health and human services task force on Alzheimer's disease. Neurology 34, 939–944. doi: 10.1212/WNL.34.7.939

PubMed Abstract | CrossRef Full Text | Google Scholar

McLachlan, D. R. C., Bergeron, C., Alexandrov, P. N., Walsh, W. J., Pogue, A. I., Percy, M. E., et al. (2019). Aluminum in neurological and neurodegenerative disease. Mol. Neurobiol. 56, 1531–1538. doi: 10.1007/s12035-018-1441-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Merlet, I., Birot, G., Salvador, R., Molaee-Ardekani, B., Mekonnen, A., Soria-Frish, A., et al. (2013). From oscillatory transcranial current stimulation to scalp EEG changes: a biophysical and physiological modeling study. PLoS ONE 8:e57330. doi: 10.1371/journal.pone.0057330

PubMed Abstract | CrossRef Full Text | Google Scholar

Meyer, M. A., Caccia, A., Martinez, D., and Mingos, M. A. (2018). PET imaging of 18F-florbetapir in cognitively impaired individuals: lack of activity within the cerebellar cortex. Neurol Int. 10:7666. doi: 10.4081/ni.2018.7666

PubMed Abstract | CrossRef Full Text | Google Scholar

Morris, E., Chalkidou, A., Hammers, A., Peacock, J., Summers, J., and Keevil, S. (2016). Diagnostic accuracy of (18)F amyloid PET tracers for the diagnosis of Alzheimer's disease: a systematic review and meta-analysis. Eur. J. Nucl. Med. Mol. Imag. 43, 374–385. doi: 10.1007/s00259-015-3228-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Moses, W. W. (2011). Fundamental limits of spatial resolution in PET. Nucl. Instrum. Methods Phys. Res. A. 648(Suppl. 1), S236–S240. doi: 10.1016/j.nima.2010.11.092

PubMed Abstract | CrossRef Full Text | Google Scholar

Murray, M. E., Lowe, V. J., Graff-Radford, N. R., Liesinger, A. M., Cannon, A., Przybelski, S. A., et al. (2015). Clinicopathologic and (11)C-Pittsburgh compound B implications of Thal amyloid phase across the Alzheimer's disease spectrum. Brain 138, 1370–1381. doi: 10.1093/brain/awv050

CrossRef Full Text | Google Scholar

Nakamura, A., Cuesta, P., Fernandez, A., Arahata, Y., Iwata, K., Kuratsubo, I., et al. (2018). Electromagnetic signatures of the preclinical and prodromal stages of Alzheimer's disease. Brain 141, 1470–1485. doi: 10.1093/brain/awy044

PubMed Abstract | CrossRef Full Text | Google Scholar

Ossenkoppele, R., Jansen, W. J., Rabinovici, G. D., Knol, D. L., van der Flier, W. M., van Berckel, B. N., et al. (2015). Prevalence of amyloid PET positivity in dementia syndromes: a meta-analysis. JAMA 313, 1939–1949. doi: 10.1001/jama.2015.4669

PubMed Abstract | CrossRef Full Text | Google Scholar

Osswald, G. (2018). BioArctic Announces Positive Topline Results of BAN2401 Phase 2b at 18 Months in Early Alzheimer's Disease. BioArctic Press Release.

Panza, F., Lozupone, M., Dibello, V., Greco, A., Daniele, A., Seripa, D., et al. (2019). Are antibodies directed against amyloid-beta (Abeta) oligomers the last call for the Abeta hypothesis of Alzheimer's disease? Immunotherapy 11, 3–6. doi: 10.2217/imt-2018-0119

PubMed Abstract | CrossRef Full Text | Google Scholar

Petkoski, S., Palva, J. M., and Jirsa, V. K. (2018). Phase-lags in large scale brain synchronization: methodological considerations and in-silico analysis. PLOS Comput. Biol. 14:e1006160. doi: 10.1371/journal.pcbi.1006160

PubMed Abstract | CrossRef Full Text | Google Scholar

Petkoski, S., Spiegler, A., Proix, T., Aram, P., Temprado, J.-J., and Jirsa, V. K. (2016). Heterogeneity of time delays determines synchronization of coupled oscillators. Phys. Rev. 94:012209. doi: 10.1103/PhysRevE.94.012209

PubMed Abstract | CrossRef Full Text | Google Scholar

Pink, J., O'Brien, J., Robinson, L., and Longson, D. (2018). Dementia: assessment, management and support: summary of updated NICE guidance. BMJ. 361:k2438. doi: 10.1136/bmj.k2438

PubMed Abstract | CrossRef Full Text | Google Scholar

Prasansuklab, A., and Tencomnao, T. (2013). Amyloidosis in Alzheimer's disease: the toxicity of Amyloid Beta (Aβ), mechanisms of its accumulation and implications of medicinal plants for therapy. Evid. Based Complement. Alternat. Med. 2013:413808. doi: 10.1155/2013/413808

CrossRef Full Text | Google Scholar

Proix, T., Bartolomei, F., Guye, M., and Jirsa, V. K. (2017). Individual brain structure and modelling predict seizure propagation. Brain 140, 641–654. doi: 10.1093/brain/awx004

PubMed Abstract | CrossRef Full Text | Google Scholar

Ren, S.-Q., Yao, W., Yan, J.-Z., Jin, C., Yin, J.-J., Yuan, J., et al. (2018). Amyloid β causes excitation/inhibition imbalance through dopamine receptor 1-dependent disruption of fast-spiking GABAergic input in anterior cingulate cortex. Sci. Rep. 8:302. doi: 10.1038/s41598-017-18729-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Reuter, M., Schmansky, N. J., Rosas, H. D., and Fischl, B. (2012). Within-subject template estimation for unbiased longitudinal image analysis. Neuroimage 61, 1402–1418. doi: 10.1016/j.neuroimage.2012.02.084

PubMed Abstract | CrossRef Full Text | Google Scholar

Rice, D. M., Buchsbaum, M. S., Starr, A., Auslander, L., Hagman, J., and Evans, W. J. (1990). Abnormal EEG slow activity in left temporal areas in senile dementia of the Alzheimer type. J. Gerontol. 45, M145–M151. doi: 10.1093/geronj/45.4.M145

PubMed Abstract | CrossRef Full Text | Google Scholar

Ripoli, C., Cocco, S., Li Puma, D. D., Piacentini, R., Mastrodonato, A., Scala, F., et al. (2014). Intracellular accumulation of amyloid-β (Aβ) protein plays a major role in Aβ-induced alterations of glutamatergic synaptic transmission and plasticity. J. Neurosci. 34, 12893–12903. doi: 10.1523/JNEUROSCI.1201-14.2014

PubMed Abstract | CrossRef Full Text | Google Scholar

Ritter, P., Schirner, M., McIntosh, A. R., and Jirsa, V. K. (2013). The virtual brain integrates computational modeling and multimodal neuroimaging. Brain Connect. 3, 121–145. doi: 10.1089/brain.2012.0120

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson, E. C., Jbabdi, S., Glasser, M. F., Andersson, J., Burgess, G. C., Harms, M. P., et al. (2014). MSM: a new flexible framework for multimodal surface matching. Neuroimage 100, 414–426. doi: 10.1016/j.neuroimage.2014.05.069

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson, J. L., Corrada, M. M., Kovacs, G. G., Dominique, M., Caswell, C., Xie, S. X., et al. (2018a). Non-Alzheimer's contributions to dementia and cognitive resilience in The 90+ Study. Acta Neuropathol. 136, 377–388. doi: 10.1007/s00401-018-1872-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson, J. L., Lee, E. B., Xie, S. X., Rennert, L., Suh, E., Bredenberg, C., et al. (2018b). Neurodegenerative disease concomitant proteinopathies are prevalent, age-related and APOE4-associated. Brain 141, 2181–2193. doi: 10.1093/brain/awy146

PubMed Abstract | CrossRef Full Text | Google Scholar

Roche (2019). Roche to Discontinue Phase III CREAD 1 and 2 Clinical Studies of Crenezumab in Early Alzheimer's Disease (AD)-Other Company Programmes in AD Continue 2019 [updated 04/01/2019]. Available online at: https://www.roche.com/dam/jcr:e3181d56-8cac-4db8-a7d4-2f883ee2847c/en/20190130-MR_CREN_EN.pdf

Sadigh-Eteghad, S., Sabermarouf, B., Majdi, A., Talebi, M., Farhoudi, M., and Mahmoudi, J. (2015). Amyloid-beta: a crucial factor in Alzheimer's disease. Med. Principles Pract. 24, 1–10. doi: 10.1159/000369101

PubMed Abstract | CrossRef Full Text | Google Scholar

Salloway, S., Honigberg, L. A., Cho, W., Ward, M., Friesenhahn, M., Brunstein, F., et al. (2018). Amyloid positron emission tomography and cerebrospinal fluid results from a crenezumab anti-amyloid-beta antibody double-blind, placebo-controlled, randomized phase II study in mild-to-moderate Alzheimer's disease (BLAZE). Alzheimer's Res. Ther. 10:96. doi: 10.1186/s13195-018-0424-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Sanz Leon, P., Knock, S. A., Woodman, M. M., Domide, L., Mersmann, J., McIntosh, A. R., et al. (2013). The Virtual Brain: a simulator of primate brain network dynamics. Front. Neuroinform. 7:10. doi: 10.3389/fninf.2013.00010

PubMed Abstract | CrossRef Full Text

Sanz-Leon, P., Knock, S. A., Spiegler, A., and Jirsa, V. K. (2015). Mathematical framework for large-scale brain network modeling in The Virtual Brain. Neuroimage 111, 385–430. doi: 10.1016/j.neuroimage.2015.01.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Schirner, M., Rothmeier, S., Jirsa, V. K., McIntosh, A. R., and Ritter, P. (2015). An automated pipeline for constructing personalized virtual brains from multimodal neuroimaging data. Neuroimage 117, 343–357. doi: 10.1016/j.neuroimage.2015.03.055

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmolesky, M. (2005). “The primary visual cortex,” in Webvision: The Organization of the Retina and Visual System [Internet], eds H. Kolb, E. Fernandez, and R. Nelson (Salt Lake City, UT: University of Utah Health Sciences Center).

Google Scholar

Schöll, M., Lockhart Samuel, N., Schonhaut Daniel, R., O'Neil James, P., Janabi, M., Ossenkoppele, R., et al. (2016). PET imaging of Tau deposition in the aging human brain. Neuron 89, 971–982. doi: 10.1016/j.neuron.2016.01.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Schüz, A., and Braitenberg, V. (ed) (2002). “The human cortical white matter: quantitative aspects of cortico-cortical long- range connectivity,” in Cortical Areas: Unity and Diversity (London; New York, NY: Taylor & Davis), 377f. doi: 10.4324/9780203219911_chapter_16

CrossRef Full Text | Google Scholar

Selkoe, D. J., and Hardy, J. (2016). The amyloid hypothesis of Alzheimer's disease at 25 years. EMBO Mol. Med. 8, 595–608. doi: 10.15252/emmm.201606210

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, R. E., Tournier, J.-D., Calamante, F., and Connelly, A. (2012). Anatomically-constrained tractography: improved diffusion MRI streamlines tractography through effective use of anatomical information. NeuroImage 62, 1924–1938. doi: 10.1016/j.neuroimage.2012.06.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, R. E., Tournier, J. D., Calamante, F., and Connelly, A. (2015). SIFT2: enabling dense quantitative assessment of brain white matter connectivity using streamlines tractography. Neuroimage 119, 338–351. doi: 10.1016/j.neuroimage.2015.06.092

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, S. M., Jenkinson, M., Woolrich, M. W., Beckmann, C. F., Behrens, T. E., Johansen-Berg, H., et al. (2004). Advances in functional and structural MR image analysis and implementation as FSL. Neuroimage 23(Suppl. 1), S208–S219. doi: 10.1016/j.neuroimage.2004.07.051

PubMed Abstract | CrossRef Full Text | Google Scholar

Song, Z., Insel, P. S., Buckley, S., Yohannes, S., Mezher, A., Simonson, A., et al. (2015). Brain amyloid-beta burden is associated with disruption of intrinsic functional connectivity within the medial temporal lobe in cognitively normal elderly. J. Neurosci. 35, 3240–3247. doi: 10.1523/JNEUROSCI.2092-14.2015

PubMed Abstract | CrossRef Full Text | Google Scholar

Sotero, R. C., and Trujillo-Barreto, N. J. (2008). Biophysical model for integrating neuronal activity, EEG, fMRI and metabolism. NeuroImage 39, 290–309. doi: 10.1016/j.neuroimage.2007.08.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Sotero, R. C., Trujillo-Barreto, N. J., Iturria-Medina, Y., Carbonell, F., and Jimenez, J. C. (2007). Realistically coupled neural mass models can generate EEG rhythms. Neural. Comput. 19, 478–512. doi: 10.1162/neco.2007.19.2.478

PubMed Abstract | CrossRef Full Text | Google Scholar

Spiegler, A. (2011). Dynamics of Biologically Informed Neural Mass Models of the Brain. Ilmenau: Universitätsverlag Ilmenau.

Google Scholar

Spiegler, A., and Jirsa, V. (2013). Systematic approximations of neural fields through networks of neural masses in the virtual brain. NeuroImage 83, 704–725. doi: 10.1016/j.neuroimage.2013.06.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Spiegler, A., Kiebel, S. J., Atay, F. M., and Knösche, T. R. (2010). Bifurcation analysis of neural mass models: impact of extrinsic inputs and dendritic time constants. NeuroImage 52, 1041–1058. doi: 10.1016/j.neuroimage.2009.12.081

PubMed Abstract | CrossRef Full Text | Google Scholar

Spiegler, A., Knösche, T. R., Schwab, K., Haueisen, J., and Atay, F. M. (2011). Modeling brain resonance phenomena using a neural mass model. PLOS Comput. Biol. 7:e1002298. doi: 10.1371/journal.pcbi.1002298

PubMed Abstract | CrossRef Full Text | Google Scholar

Stefanescu, R. A., and Jirsa, V. K. (2008). A low dimensional description of globally coupled heterogeneous neural networks of excitatory and inhibitory neurons. PLoS Comput. Biol. 4:e1000219. doi: 10.1371/journal.pcbi.1000219

PubMed Abstract | CrossRef Full Text | Google Scholar

Stefanovski, L., Ghani, A., McIntosh, A. R., and Ritter, P. (2016). Linking connectomics and dynamics in the human brain. e-Neuroforum 7, 64–70. doi: 10.1007/s13295-016-0027-1

CrossRef Full Text | Google Scholar

Storck, S. E., and Pietrzik Claus, U. (2018). The Blood brain-barrier and its role in Alzheimer's disease. Neuroforum 24:A197. doi: 10.1515/nf-2018-A014

CrossRef Full Text | Google Scholar

Strogatz, S. H. (2015). Nonlinear Dynamics and Chaos: With Applications to Physics, Biology, Chemistry, and Engineering: Second Edition. Boulder, CO: Westview Press, a member of the Perseus Books Group.

Google Scholar

Tadel, F., Baillet, S., Mosher, J. C., Pantazis, D., and Leahy, R. M. (2011). Brainstorm: a user-friendly application for MEG/EEG analysis. Comput. Intell. Neurosci. 2011:13. doi: 10.1155/2011/879716

PubMed Abstract | CrossRef Full Text | Google Scholar

Takatori, S., Wang, W., Iguchi, A., and Tomita, T. (2019). Genetic risk factors for alzheimer disease: emerging roles of microglia in disease pathomechanisms. Adv. Exp. Med. Biol. 1118, 83–116. doi: 10.1007/978-3-030-05542-4_5

PubMed Abstract | CrossRef Full Text | Google Scholar

Tapia-Rojas, C., Cabezas-Opazo, F., Deaton, C. A., Vergara, E. H., Johnson, G. V. W., and Quintanilla, R. A. (2019). It's all about tau. Prog. Neurobiol. 175, 54–76. doi: 10.1016/j.pneurobio.2018.12.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Thal, D. R., Rub, U., Orantes, M., and Braak, H. (2002). Phases of a beta-deposition in the human brain and its relevance for the development of AD. Neurology 58, 1791–1800. doi: 10.1212/WNL.58.12.1791

PubMed Abstract | CrossRef Full Text | Google Scholar

Thomas, B. A., Cuplov, V., Bousse, A., Mendes, A., Thielemans, K., Hutton, B. F., et al. (2016). PETPVC: a toolbox for performing partial volume correction techniques in positron emission tomography. Phys. Med. Biol. 61, 7975–7993. doi: 10.1088/0031-9155/61/22/7975

PubMed Abstract | CrossRef Full Text | Google Scholar

Tournier, J. D., Calamante, F., and Connelly, A. (2007). Robust determination of the fibre orientation distribution in diffusion MRI: non-negativity constrained super-resolved spherical deconvolution. Neuroimage 35, 1459–1472. doi: 10.1016/j.neuroimage.2007.02.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Tournier, J. D., Calamante, F., and Connelly, A. (2010). Improved probabilistic streamlines tractography by 2nd order integration over fibre orientation distributions. Proc. Intl. Soc. Mag. Reson. Med. 18:1670.

Google Scholar

Tournier, J. D., Calamante, F., and Connelly, A. (2013). Determination of the appropriate b value and number of gradient directions for high-angular-resolution diffusion-weighted imaging. NMR Biomed. 26, 1775–1786. doi: 10.1002/nbm.3017

PubMed Abstract | CrossRef Full Text | Google Scholar

Ulrich, D. (2015). Amyloid-beta impairs synaptic inhibition via GABA(A) receptor endocytosis. J. Neurosci. 35, 9205–9210. doi: 10.1523/JNEUROSCI.0950-15.2015

PubMed Abstract | CrossRef Full Text | Google Scholar

Veraart, J., Novikov, D. S., Christiaens, D., Ades-aron, B., Sijbers, J., and Fieremans, E. (2016). Denoising of diffusion MRI using random matrix theory. NeuroImage 142, 394–406. doi: 10.1016/j.neuroimage.2016.08.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Verret, L., Mann Edward, O., Hang Giao, B., Barth Albert, M. I., Cobos, I., Ho, K., et al. (2012). Inhibitory interneuron deficit links altered network activity and cognitive dysfunction in Alzheimer model. Cell 149, 708–721. doi: 10.1016/j.cell.2012.02.046

PubMed Abstract | CrossRef Full Text | Google Scholar

Villemagne, V. L., Ataka, S., Mizuno, T., Brooks, W. S., Wada, Y., Kondo, M., et al. (2009). High striatal amyloid beta-peptide deposition across different autosomal Alzheimer disease mutation types. Arch. Neurol. 66, 1537–1544. doi: 10.1001/archneurol.2009.285

PubMed Abstract | CrossRef Full Text | Google Scholar

Walsh, D. M., and Selkoe, D. J. (2007). A beta oligomers - a decade of discovery. J. Neurochem. 101, 1172–1184. doi: 10.1111/j.1471-4159.2006.04426.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, S., and Colonna, M. (2019). Microglia in Alzheimer's disease: a target for immunotherapy. J. Leukoc. Biol. 106, 219–227. doi: 10.1002/jlb.mr0818-319r

PubMed Abstract | CrossRef Full Text | Google Scholar

Wendling, F., Bartolomei, F., Bellanger, J. J., and Chauvel, P. (2002). Epileptic fast activity can be explained by a model of impaired GABAergic dendritic inhibition. Eur. J. Neurosci. 15, 1499–1508. doi: 10.1046/j.1460-9568.2002.01985.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wendling, F., Bellanger, J. J., Bartolomei, F., and Chauvel, P. (2000). Relevance of nonlinear lumped-parameter models in the analysis of depth-EEG epileptic signals. Biol. Cybernet. 83, 367–378. doi: 10.1007/s004220000160

PubMed Abstract | CrossRef Full Text | Google Scholar

WHO (2011). Global Health and Ageing. Aging UNIo, WHO.

Wilson, H. R., and Cowan, J. D. (1972). Excitatory and inhibitory interactions in localized populations of model neurons. Biophys. J. 12, 1–24. doi: 10.1016/S0006-3495(72)86068-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Wimo, A., Guerchet, M., Ali, G. C., Wu, Y. T., Prina, A. M., Winblad, B., et al. (2017). The worldwide costs of dementia 2015 and comparisons with 2010. Alzheimers Dement. 13, 1–7. doi: 10.1016/j.jalz.2016.07.150

PubMed Abstract | CrossRef Full Text | Google Scholar

Wimo, A., Jonsson, L., Gustavsson, A., McDaid, D., Ersek, K., Georges, J., et al. (2011). The economic impact of dementia in Europe in 2008-cost estimates from the Eurocode project. Int. J. Geriatr. Psychiatry 26, 825–832. doi: 10.1002/gps.2610

PubMed Abstract | CrossRef Full Text | Google Scholar

Wong, K.-F., and Wang, X.-J. (2006). A recurrent network mechanism of time integration in perceptual decisions. J. Neurosci. 26, 1314–1328. doi: 10.1523/JNEUROSCI.3733-05.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Woolrich, M. W., Jbabdi, S., Patenaude, B., Chappell, M., Makni, S., Behrens, T., et al. (2009). Bayesian analysis of neuroimaging data in FSL. Neuroimage 45(1 Suppl.), S173–S186. doi: 10.1016/j.neuroimage.2008.10.055

PubMed Abstract | CrossRef Full Text | Google Scholar

Wright, J. J., and Liley, D. T. J. (2010). Dynamics of the brain at global and microscopic scales: neural networks and the EEG. Behav. Brain Sci. 19, 285–295. doi: 10.1017/S0140525X00042679

CrossRef Full Text | Google Scholar

Xu, J., Wang, J., Wimo, A., Fratiglioni, L., and Qiu, C. (2017). The economic burden of dementia in China, 1990-2030: implications for health policy. Bull. World Health Organ. 95, 18–26. doi: 10.2471/BLT.15.167726

PubMed Abstract | CrossRef Full Text | Google Scholar

Zetterberg, L. H., Kristiansson, L., and Mossberg, K. (1978). Performance of a model for a local neuron population. Biol. Cybernet. 31, 15–26. doi: 10.1007/BF00337367

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, J., Yu, W., Zhang, M., Tian, X., Li, Y., and Lu, Y. (2019). Imbalance of microglial TLR4/TREM2 in LPS-treated APP/PS1 transgenic mice: a potential link between Alzheimer's disease and systemic inflammation. Neurochem. Res. 44, 1138–1151. doi: 10.1007/s11064-019-02748-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Zimmermann, J., Perry, A., Breakspear, M., Schirner, M., Sachdev, P., Wen, W., et al. (2018). Differentiation of Alzheimer's disease based on local and global parameters in personalized Virtual Brain models. NeuroImage Clin. 19, 240–251. doi: 10.1016/j.nicl.2018.04.017

PubMed Abstract | CrossRef Full Text | Google Scholar