Computational models in Electroencephalography
Katharina Glomb, Joana Cabral, Anna Cattani, Alberto Mazzoni, Ashish Raj, Benedetta Franceschiello
NNoname manuscript No. (will be inserted by the editor)
Computational models in Electroencephalography
Katharina Glomb · Joana Cabral · Anna Cattani · Alberto Mazzoni · Ashish Raj · Benedetta Franceschiello
Received: date / Accepted: date
Abstract
Computational models lie at the intersectionof basic neuroscience and healthcare applications be-cause they allow researchers to test hypotheses in silico
K. GlombDepartment of Radiology,Lausanne University Hospital and University of Lausanne(CHUV-UNIL),Lausanne, SwitzerlandE-mail: [email protected]. CabralLife and Health Sciences Research Institute (ICVS),University of Minho,PortugalA. CattaniDepartment of Psychiatry, University of Wisconsin,Wisconsin, USADepartment of Biomedical and Clinical Sciences ’Luigi Sacco’,University of Milan,Milan, ItalyA. MazzoniThe BioRobotics Institute,Scuola Superiore Sant’Anna,Pisa, ItalyA. RajSchool of Medicine, UCSF,San Francisco, USB. FranceschielloDepartment of Ophthalmology,Hopital Ophthalmic Jules Gonin, FAA,Lausanne, SwitzerlandLaboratory for Investigative Neurophysiology,Department of Radiology,Lausanne University Hospital and University of Lausanne(CHUV-UNIL),Lausanne, SwitzerlandE-mail: [email protected] and predict the outcome of experiments and interactionsthat are very hard to test in reality. Yet, what is meantby “computational model” is understood in many differ-ent ways by researchers in different fields of neuroscienceand psychology, hindering communication and collabo-ration. In this review, we point out the state of the artof computational modeling in Electroencephalography(EEG) and outline how these models can be used tointegrate findings from electrophysiology, network-levelmodels, and behavior. On the one hand, computationalmodels serve to investigate the mechanisms that gen-erate brain activity, for example measured with EEG,such as the transient emergence of oscillations at dif-ferent frequency bands and/or with different spatialtopographies. On the other hand, computational modelsserve to design experiments and test hypotheses in silico .The final purpose of computational models of EEG isto obtain a comprehensive understanding of the mech-anisms that underlie the EEG signal. This is crucialfor an accurate interpretation of EEG measurementsthat may ultimately serve in the development of novelclinical applications.
Keywords: electroencephalography, computationalmodeling, multiscale modeling, clinical applications
Electroencephalography (EEG) has applications in manyfields, spanning from basic neuroscientific research toclinical domains. However, despite the technological ad-vances in recording precision, the full potential of EEGis currently not being exploited. One possible way to doso is to use computational models in order to integratefindings from electrophysiology, network-level models(the level of neuroimaging), and behavior (Franceschiello a r X i v : . [ q - b i o . N C ] S e p t al., 2018, 2019). A model is defined in terms of a setof equations which describe the relationships betweenvariables. Importantly, models exist for different spatialscales (Varela et al., 2001; Deco et al., 2008a), span-ning from the single cell spike train up to macroscopicoscillations. The equations are used to simulate howeach variable changes over time, or, in rare cases, tofind analytical solutions for the relationships among thevariables. The dynamics of the resulting time series arealso influenced by a set of parameters, which can eitherbe estimated from available data - for example, a modelwhich simulates the firing of a certain neuron type couldcontain a time constant estimated from recordings onthat type of neuron in rodents - or its value can be variedsystematically in an exploratory manner. The goal is toproduce time series of variables that can be compared toreal data. In particular, one can simulate perturbationsto brain activity, be it sensory stimulation, a therapeuticintervention like DBS or a drug, or a structural changedue to the onset of a pathology, like neurodegenerationor a lesion, and predict what would be the resultingalterations observed in neural and clinical data.An important application of EEG models is in theclinical domain. Psychiatric and neurological disordersimpact a growing portion of the population, both as pa-tients and caregivers, and with an enormous cost - botheconomical and humanitarian - to healthcare systemsworldwide (Steel et al., 2014; Vigo et al., 2016; Feiginet al., 2019). One of the main obstacles in advancingpatient care is the lack of individual diagnosis, prognosis,and treatment planning (Wium-Andersen et al., 2017).Computational models can be adapted to the individualby setting their parameters according to available data(i.e. either setting the parameter directly, if it is measur-able, or looking for the parameter value which results intime series whose dynamics match recorded data). Theadjusted parameter(s) can then be related to clinicalmarkers, symptoms, and behavior, allowing for exampleto discriminate between pathologies. Using models inthis personalized manner could provide additional di-agnostic features in the form of model parameters andmodel output, eventually assisting clinicians in diagnosisand treatment planning.Another obstacle is a general lack of scientific knowl-edge of disease mechanisms, including the mechanismsby which therapies exert their effect. As an example,deep brain stimulation (DBS) is a highly effective treat-ment for advanced Parkinsonism, in which electric pulsesare delivered directly to certain deep brain structuresvia permanently implanted electrodes. Yet, it is largelyunknown how exactly the applied stimulation managesto suppress motor symptoms such as tremor (Chikenand Nambu, 2016). This is also because the way in which motor symptoms result from the degradation ofdopaminergic neurons in the substantia nigra is notfully understood (McGregor and Nelson, 2019). Besidesanimal models - which have their own ethical issues - insilico models are an indispensable tool for understand-ing brain disorders. Combining data available from apatient or group of patients with knowledge and hy-potheses about mechanisms, a model can be generatedwhich can help test these hypotheses.Last but not least, models are much cheaper thananimal testing or clinical trials. While models will notreplace these approaches - at least not in the foresee-able future - they could help to formulate more specifichypotheses and thus, lead to smaller-scale experiments.Collecting invasive data is not generally possible in hu-mans. EEG is an extremely versatile technology whichallows non-invasive recording of neural activity in behav-ing humans. EEG is a cheap and portable technology,particularly compared to (f)MRI and MEG. Apart fromthese cost-efficiency considerations, EEG, like MEG, isa direct measure of the electromagnetic fields generatedby the brain, and allows millisecond-precision record-ings, thus giving access to rich aspects of brain functionwhich can inform models in a way that e.g. fMRI cannot(see section 2 for more details). In general, using dif-ferent complementary sources of data to construct andvalidate a model will lead to better model predictions,as each recording technique has its own strengths andweaknesses, and a multimodal approach can balancethem.In our opinion, there are mostly two reasons whyEEG has not been used more extensively in modelingstudies, and particularly in a clinical context. First,there are numerous technical problems which make theprocessing and interpretation of EEG data challenging.EEG - like MEG - is measured on the scalp, and theproblem of projecting this 2D-space into the 3D-brainspace arises (Michel and Brunet, 2019). While multiplesolutions exist for this inverse problem, it is unclearwhich one is the best and under which circumstances(Hassan et al., 2014; Mahjoory et al., 2017; Hedrichet al., 2017). EEG data require extensive preprocess-ing, e.g. removal of artifacts due to movements, eyeblinks, etc., but these steps are far from being stan-dardized, and many options exist. The recently startedEEG-BIDS effort (Pernet et al., 2019) is a step towardsthe direction of standardization of EEG data and shouldfacilitate, alongside with the much larger amount of pub-licly available data, studies that systematically evaluatethe impact of preprocessing steps and compare sourcereconstruction algorithms. As the interest in EEG rises,the need to resolve these issues will trigger larger effortsthat will benefit the entire community.2he second obstacle to a more routine usage of com-putational models in EEG research, which we hope toaddress in this review, is that such models usually re-quire an understanding of the mathematics involved, ifonly to be able to choose the model that is useful forthe desired application. Both variables and parametersare not always clearly related to quantities which canbe measured in a clinical or experimental context, andmore generally, models need to be set up in such a waythat they meet existing clinical demands or researchquestions.The contribution of this paper is threefold. First,this article summarizes computational approaches atdifferent spatial scales in EEG, targeting non-expertsreaders. To the best of our knowledge, this paper rep-resents the first review on this topic. Second, we willpoint out several ways in which computational modelsintegrate EEG recordings, by using biologically relevantvariables. Third, we discuss the clinical applicationsof computational models in EEG which have been de-veloped. The field is greatly expanding and containspromising advancements both from research and clinicalstandpoints. We believe that this overview will make thefield accessible for a broad audience, and indicate thenext steps required to push modeling of EEG forward. EEG is a non-invasive neuroimaging technique that mea-sures the electrical activity of the brain (Biasiucci et al.,2019). EEG recordings have been a driver of researchand clinical applications in neuroscience and neurologyfor nearly a century. EEG relies on the placement of elec-trodes on the person’s scalp, measuring the postsynapticpotentials of pyramidal neurons (Tivadar and Murray,2019; da Silva, 2013). EEG does not directly measurethe action potentials of neurons, though there are someindications of high-frequency oscillations being linked tospiking activity (Telenczuk et al., 2011). The neurotrans-mitter release generated by action potentials, whetherexcitatory or inhibitory, results in local currents at theapical dendrites that in turn lead to current sources andsinks in the extracellular space around the dendriticarbor (i.e. postsynaptic potentials, see Figure 1, bot-tom right block). EEG is generated by the local fieldpotential (LFP), a signal that reflects summed synapticactivity of local populations of neurons. In the neocortex,pyramidal neurons are generally organized perpendicu-larly to the cortical surface, with apical dendrites towardthe pial surface and axons pointing inferiorly towardsthe grey-white matter border. This alignment leads tothe electrical fields of many neurons being summed upto generate a signal that is measurable at the scalp (Tivadar and Murray, 2019). Importantly, individualneurons of these populations need to be (nearly) syn-chronously active to be detectable by EEG.As mentioned above, the electrical activity of thebrain is recorded by means of electrodes, made of con-ductive materials, placed at the scalp. The propagationof electrical fields takes place due to the conductiveproperties of brain and head tissues, a phenomenonknown as volume conduction (Kajikawa and Schroeder,2011). The electrodes are connected to an amplifierwhich boosts the signal. Due to the biophysical natureof what is measured, i.e. a voltage - the difference ofpotential able to move charges from one site to the other- EEG records the differential measurements betweenan electrode at a specific position on the scalp and areference site. Common analyses in EEG are the studyof local phenomena such as peaks at specific latencies orscalp sites (event-related potential, ERPs); or the studyof topography, i.e. the shape of the electric field at thescalp, which represents a global brain signature (Mur-ray et al., 2008). EEG is known for its high temporalresolution. The biggest pitfalls of the technique are, onthe other hand, the low spatial resolution and signal tonoise ratio. A clear and exhaustive walk through thesetopics as well as an overview of strengths and pitfalls ofusing EEG is contained in Tivadar and Murray (2019)and for non-experts of the field in Biasiucci et al. (2019).Despite being a measurement of the scalp activity,EEG can reveal the underlying neurophysiological mech-anisms of the brain, and that is what classifies it asbrain imaging tool. The estimation of the loci of activesources for the recorded brain activity at the scalp iscalled source reconstruction (Michel et al., 2004). How-ever, the loci can belong to areas not necessarily belowthe considered electrode, a pitfall caused by volumeconduction. Source reconstruction is a mathematicalill-posed inverse problem, as the solution is not unique.However, the addition of biophysical constraints to theinverse problem allows to retrieve a solution, whichhas been validated by means of intracranial recordings(Michel and Murray, 2012). Having obtained the sourceactivity, one can estimate the functional connectivitybetween the sources, i.e. the statistical dependenciesbetween brain areas, assumed to indicate their interac-tions (see also table 1). This can then be complementedwith neuroanatomical/structural connectivity (table 1),which estimates white matter connections between brainareas.Computational models stand at the interface be-tween the physiology of neurons at different scales (singleneuron, population, macro-scale) and perceptual behav-ior. EEG would greatly benefit from the integrationof in-silico simulations, as computational models could3omplement both the neurophysiological and behavioralinterpretations of EEG recordings. In the following sec-tions, we will discuss different types of computationalmodels, i.e. the different scales at which the neural ac-tivity is simulated, how such models can be integratedin the analysis of EEG signals, and how such modelshave been used in new clinical applications.
A straightforward classification of computational modelsfor EEG can be done based on the different scales ofneurophysiological activity they integrate. For instance,we can distinguish three types of models (Figure 2 A):1. microscopic models on the level of single cells andmicro-circuits;2. mesoscopic models on the level of neural masses andneural fields;3. macroscopic models taking into account the connec-tome/white matter.The integration of computational models has greatlyadvanced the field of applications of EEG, both forresearch and clinical purposes.3.1 Computational models for EEG on the level ofsingle cells and microcircuitsThe purpose of this level of modeling is to addressthe origin of the EEG signal by investigating the rela-tionship between its features and electrophysiologicalmechanisms (Figure 1, column A) with the tools of com-putational neuroscience. As detailed above, the EEGsignal recorded from the scalp is the result of the spatialintegration of the potential fluctuations in the extracel-lular medium. The EEG signal is mainly caused by thelocal field potential (LFP), while LFP is mainly drivenby synaptic activity (Logothetis, 2003; Buzski et al.,2012) and volume conduction (Kajikawa and Schroeder,2011). From the experimental standpoint, local networkactivity is usually measured as LFP (mainly in vivo -and rarely in vitro - animal data). By virtue of superpo-sition, fluctuations in the LFP, and EEG more generally,are signatures of correlated neural activity (Pesaranet al., 2018). Cellular and microcircuit modeling arethus aimed at understanding the neurophysiological un-derpinnings of these correlations and the role played bycell types, connectivity and other properties in shapingthe collective activity of neurons.A primary goal of EEG modeling at the microscopicscale is on the one hand to predict the EEG signal generated by the summation of local dynamics on themicroscopic scale and, on the other hand, to reconstructthe microscopic neural activity underlying the observedEEG. The first goal is far from being achieved, and thesecond is ill-posed due to the number of possible circuitand cellular combinations at the source level leading tosimilar EEGs. Implicit to these goals is to understandhow features of neural circuits, such as the architecture,synapses and cell types, contribute to the generation ofelectromagnetic fields and their properties in a bottom-up fashion. Despite key insights, many shortcomingslimit the interpretability of microcircuit models andthe establishment of a one-to-one correspondence withEEG data. For instance, the contribution of spikingactivity and correlated cellular fluctuations to LFPs andEEG power spectra remains unclear. Most microcircuitmodels characterize the net local network activity - usedas a proxy for EEG - using the average firing rate or viathe mean somatic membrane potential taken amongstpopulations of cells (of various types). Other studieshave used a heuristic approach and approximated theEEG signal as a linear combination of somatic membranepotentials with random coefficients to account for bothconduction effects and observational noise (Herrmannet al., 2016; Lefebvre et al., 2017). As such, microcircuitmodel predictions and experimental data cannot alwaysbe compared directly.Cellular multicompartmental models, which often-times take cellular morphology and spatial configurationinto consideration, are based on the celebrated Hodgkin-Huxley equations, which describe the temporal evolutionof ionic flux across neuronal membranes (see Catterallet al. (2012), for a recent review). Such conductance-based models, which possess explicit and spatially dis-tributed representations for cellular potentials, facilitatethe prediction and/or comparison with LFP recordings.In contrast, single compartment models are difficultto interpret: while more abstract single compartmentmodels such as Poisson neurons or integrate-and-firemodels (Figure 2 A, left) are often used for their relativetractability and computational efficiency to constructmore elaborate microcircuit models, they generally lackthe neurophysiological richness to estimate EEG traces.Despite this, several computational advancements inrecent years investigated how networks of integrate-and-fire neurons generate LFPs, clarifying the microscopicdynamics reflected in the EEG signal (Mazzoni et al.,2015, 2008, 2010, 2011; Deco et al., 2008b; Buehlmannand Deco, 2008; Barbieri et al., 2014). Such approacheshave been used to understand the formation of correlatedactivity patterns in the hippocampus (e.g. oscillations),and their associated spectral fingerprints in the LFP(Chatzikalymniou and Skinner, 2018). Furthermore, a4 ig. 1
Electrophysiology of neural activity and EEG at different scales. ig. 2 A - Illustration of computational models at the three scales treated here. Microscopic scale : Simple example of two( i = 1 ,
2) leaky integrate-and-fire (LIF) neurons coupled together, a pyramidal neuron making an excitatory synapse to theinterneuron, which in turn makes an inhibitory synapse to the pyramidal cell. This minimal circuit implements feedbackinhibition, as the pyramidal cell, when activated, will excite the interneuron, which in turn will inhibit it. In the equation, V i is the membrane potential of each of the two cells i = 1 , V L is the leak, or resting potential of the cells; R is a constantcorresponding to the membrane resistance; I i is the synaptic input that each cell receives from the other, and possibly backgroundinput; τ is the time constant determining how quickly V i decays. The model is simulated by setting a firing threshold, at which,when reached, a spike is recorded and V i is reset to V L . Mesoscopic scale : The Wilson-Cowan-model, in which an excitatory( E ) and an inhibitory ( I ) population are coupled together. The mean field equations describe the mean activity of a largenumber of neurons. f E and f I are sigmoid transfer functions whose values indicate how many neurons in the population reachfiring threshold, and h E / h I are external inputs like background noise. w EE and w II are constants correponding to the strengthof self-excitation/inhibition, and w EI and w IE the strength of synaptic coupling between populations. Macroscopic scale : Inorder to simulate long-range interactions between cortical and even subcortical areas, brain network models couple togethermany mesoscopic (“local”) models using the connection weights defined in the empirical structural connectivity matrix C . Theexample equation defines the Kuramoto model, in which the phase θ n of each node n is used as a summary of its oscillatoryactivity around its natural frequency ω . Each node’s phase depends on the phases of connected nodes p taking into accountthe time delay τ np , defined by the distances between nodes n and p . k is a global scaling parameter controlling the strengthof internode connections. B - Illustration of a typical modeling approach at the macroscopic scale. Activity is simulated foreach node using the defined macroscopic model, e.g. the Kuramoto model from panel A, right. The feature of interest is thencomputed from this activity. Shown here is the functional connectivity, e.g., phase locking values between nodes (table 1). Thiscan then be compared to the empirical functional connectivity matrix computed in exactly the same way from experimentaldata, e.g. by correlating the entries of the matrix. The model fit can be determined depending on parameters of the model, e.g.the scaling parameter k or the unit speed, here indicated with “tau”. w IE and w EI in Figure 2 A, middle), both within andbetween populations. Also it is possible to modify timeconstants (which govern the decay of activity in thelocal populations) or the strength of background noise.Changing these parameters in silico can be interpretedbiologically. For example, in Bojak and Liley (2005),the authors describe how a modified neural field modelreproduces EEG spectra recorded during anaesthesia.The strength of inputs from the thalamus to the corti-cal neural populations was varied within a biologicallyplausible range.Neural mass and neural field models are able to repro-duce a range of dynamical behaviors that are observed inEEG, like oscillations in typical EEG frequency bands(David and Friston, 2003), phase-amplitude-coupling(Onslow et al., 2014; Sotero, 2016), evoked responses(Jansen et al., 1993; Jansen and Rit, 1995; David et al.,2005), and allows to model the EEG spectrum (Davidand Friston, 2003; Bojak and Liley, 2005; Moran et al.,2007).By coupling together more than one model/set ofpopulations, one can start investigating the effect thatdelays have on neural activity (Jirsa and Haken, 1996).In fact, Jansen & Rit (Jansen and Rit, 1995) coupledtogether two neural mass models in order to simulatethe effect of interactions between cortical columns ontheir activity.Often, activity simulated by mean field models is as-sumed to be related to local field potentials (Liley et al.,2002). However, models are usually set up such thatthe local field potential derives directly from the meanfiring rate. In this way, an important aspect that under-lies the EEG signal is neglected, namely, the synchrony(coherence) of the firing within a neural population (asopposed to synchrony between populations, which canbe studied using e.g. instantaneous phase differences(Breakspear et al., 2004)). Phenomena such as event-related synchronization and -desynchronization resultfrom a change in this synchrony rather than from achange in firing rate. Recent models (Byrne et al., 2017,2020) propose therefore a link between the firing rateand the Kuramoto order parameter, which is a measureof how dispersed firing is within a population. 3.3 Macroscopic computational models for EEG takinginto account the connectomeIn this section, we review existing literature on macro-scopic computational models that take into account theconnectome and discuss their potential to reveal thegenerative mechanisms of the macroscopic brain activitypatterns detected with EEG and MEG (Figure 1, col-umn C). We will use the term “brain network models”(BNM) in order to clearly distinguish this frameworkfrom other approaches to whole-brain modeling (Break-spear, 2017), e.g. using neural field models (Jirsa andHaken, 1996; Robinson et al., 1997; Coombes et al.,2007) or expansions of the thalamocortical models dis-cussed above (Robinson et al., 2001; Freyer et al., 2011).We will also leave aside the large body of literature ondynamic causal modeling (DCM) (Kiebel et al., 2008;Pinotsis et al., 2012), as this deserves a more detailedreview than the scope of this paper can provide. Brain network models.
In recent years, the interest inthe human connectome has experienced a boom, creat-ing the prolific and successful field of “connectomics”.In the framework of connectomics, the brain is concep-tualized as a network made up of nodes and edges. Eachnode represents a brain region, and nodes are coupledtogether according to a weighted matrix representing thewiring structure of the brain (Figure 2 A, right). Thisso-called structural connectivity matrix (SC) is derivedfrom white matter fiber bundles which connect distantbrain regions (Behrens et al., 2003; Zhang et al., 2010;Hagmann et al., 2008; Sepulcre et al., 2010; Wedeenet al., 2012) and are measured using diffusion weightedmagnetic resonance imaging (dMRI) (table 1). The setof all fiber bundles is called the connectome (Sporns,2011). By coupling brain regions together according tothe weights in the SC, the activity generated in eachregion depends also on the activity propagated fromother regions along the connections given by the SC.BNMs are used to study the role of structural con-nectivity in shaping brain activity patterns. Because thisis a complex problem that involves the entire brain, it isimportant to find a balance between realism and reduc-tion, so that useful predictions can be made.In practice,a common simplification is to assume that all brain re-gions are largely identical in their dynamical properties(Passingham et al., 2002). This reductionist approachkeeps the number of parameters at a manageable leveland still allows to investigate how collective phenomenaemerge from the realistic connectivity between nodes.In other words, BNMs do not necessarily aim at maxi-mizing the fit to the empirically recorded brain signals.Rather, the goal is to reproduce specific temporal, spa-8ial or spectral features of the empirical data emergingat the macroscopic scale whose underlying mechanismsremain unclear (Figure 2 B).
Choice of local model.
In mathematical terms, brainactivity is simulated according to a system of coupleddifferential equations. The activity of each node is de-scribed by a mean-field model, such as the ones describedin section 3.2, and coupling between the mean field mod-els is parametrized by the empirical SC (Figure 1 A,right).Importantly, the type of mean-field model used atthe local level must be selected according to the hy-pothesis being tested. For example, BNMs have provedto be a powerful tool to elucidate the non-linear linkbetween the brain’s structural wiring and the functionalpatterns of brain activity captured with resting-statefunctional magnetic resonance imaging (rsfMRI) (Decoet al., 2014a, 2013a; Honey et al., 2009; Deco et al.,2009; Cabral et al., 2011). However, oscillations in fre-quency ranges important for M/EEG (2-100 Hz) areoften neglected in studies aiming at reproducing corre-lated fluctuations on the slow time scale of the fMRIsignal. Thus, despite the insights gained by BNMs tounderstand rsfMRI signal dynamics, the same modelsdo not necessarily serve to understand M/EEG signalsand vice-versa (Cabral et al., 2017).In Cabral et al. (2014), the local model employedincludes a mechanism for the generation of collectiveoscillatory signals in order to address oscillatory com-ponents of M/EEG. To model brain-wide interactionsbetween local nodes oscillating around a given naturalfrequency (in this case, 40 Hz, in the gamma frequencyrange), the Kuramoto model (Kuramoto, 2003; Yeungand Strogatz, 1999), was extended to incorporate realis-tic brain connectivity (SC) and time delays (determinedby the lengths of the fibers in the SC (see also Fingeret al. (2016); Figure 2 A, right). This model shows how,for a specific range of parameters, groups of nodes (com-munities) can temporarily synchronize at community-specific lower frequencies, obeying to universal rulesthat govern the behaviour of coupled oscillators withtime delays. Thus, the model proposed a mechanismthat explains how slow global rhythms in the alpha-and beta-range emerge from interactions of fast local(gamma) oscillations generated by neuronal networks.In contrast, Deco et al. (2019b) used a mean fieldmodel (Wilson and Cowan, 1973; Brunel and Wang,2001; Deco and Jirsa, 2012; Deco et al., 2014b), whichwas tuned not to exhibit intrinsic oscillations. Becausethe brain could thus be considered as being in a noisy,low-activity state, the number of parameters was suffi-ciently reduced to investigate how activation patterns change over time on different time scales. Time scalesincluding that of M/EEG (ten to several 100 ms) aswell as that typical for fMRI (1-3 seconds) were con-sidered, and the question was asked whether there isa time scale at which brain dynamics are particularlyrich. The authors found that the best frequency reso-lution was on the scale of 200 ms, as both the numberof co-activation patterns as well as their dynamics wererichest, compared to other resolutions.
Emerging class of harmonics-based models.
Althoughboth the described BNMs as well as DCM (dynamiccausal modeling) have a long history of success in model-ing brain activity patterns, they have high-dimensionality,and usually require local oscillators governed by region-specific or spatially-varying model parameters. Whilethis imbues such models with rich features capable ofrecreating complex behavior, they are challenging forsome clinical applications where a small set of globalfeatures might be desired to assess the effect of diseaseon network activity. Therefore recently some laborato-ries have focused on low-dimensional processes involvingdiffusion or random walks (table 1) on the structuralgraph (table 1) instead of mean-field models, providinga simpler means of simulating functional connectivity(FC). These simpler models were able to match or exceedthe predictive power of complex neural mass models orDCMs in predicting empirical FC (Abdelnour et al.,2014). Higher-order walks on graphs have also beenquite successful; typically these methods involve a seriesexpansion of the graph adjacency or Laplacian matrices(Meier et al., 2016; Becker et al., 2018) (table 1). Notsurprisingly, the diffusion and series expansion methodsare closely related, and most of these approaches maybe interpreted as special cases of each other, as demon-strated elegantly in recent studies (Robinson et al., 2016;Deslauriers-Gauthier et al., 2020; Tewarie et al., 2020).Whether using graph diffusion or series expansion,these models of spread naturally employ the so-calledeigenmodes, or harmonics, of graph adjacency or Lapla-cian matrix. Hence these methods were generalized toyield spectral graph models whereby e.g. Laplacian har-monics were sufficient to reproduce empirical FC, usingonly a few eigenmodes (Atasoy et al., 2016; Abdelnouret al., 2018). The Laplacian matrix in particular has along history in graph modeling, and its eigenmodes arethe orthonormal basis of the network and can thus repre-sent arbitrary patterns on the network (Stewart, 1999).Such spectral graph models are computationally attrac-tive due to low-dimensionality and more interpretableanalytical solutions. The SC’s Laplacian eigenmodesmay be thought of as the substrate on which functionalpatterns of the brain are established via a process of net-9ork transmission (Abdelnour et al., 2018; Atasoy et al.,2016; Robinson et al., 2016; Preti and Van De Ville,2019; Glomb et al., 2020). These models were strikinglysuccessful in replicating canonical functional networks,which are stable large scale circuits made up of function-ally distinct brain regions distributed across the cortexthat were extracted by clustering a large fMRI dataset(Yeo et al., 2011).While spectral graph models have demonstrated abil-ity to capture essential steady-state, stationary char-acteristics of real brain activity, they are limited tomodeling passive spread without oscillatory behavior.Hence they may not suitably accommodate a largerrepertoire of dynamically-varying microstates or richpower spectra at higher frequencies typically observedon EEG or MEG. Capturing this rich repertoire wouldrequire a full accounting of axonal propagation delays aswell as local neural population dynamics within graphmodels, as previously advocated (Cabral et al., 2011).Band-specific MEG resting-state networks were success-fully modeled with a combination of delayed neuralmass models and eigenmodes of the structural network(Tewarie et al., 2019), suggesting delayed interactionsin a brain’s network give rise to functional patternsconstrained by structural eigenmodes. Recently anothereffort was undertaken to characterise wide-band brainactivity using graph harmonics in closed form (i.e. re-quiring no time-domain simulations), a rarity in the fieldof computational neuroscience (Raj et al., 2020). This“spectral graph model” of brain activity produced realis-tic power spectra that could successfully predict boththe spatial as well as temporal properties of MEG empir-ical recordings (Raj et al., 2020). Intriguingly, the modelhas very few (six) parameters, all of which are globaland not dependent on local oscillations. This methodtherefore exemplifies the power of graph methods inreproducing more complex and rich repertoire of brainactivity, while keeping to a parsimonious approach thatdoes not require the kinds of high-dimensional and non-linear oscillatory models that have traditionally heldsway.
Network oscillations, captured through EEG, are thoughtto be relevant for brain functions, such as cognition,memory, perception and consciousness (Ward, 2003).Local brain regions produce oscillatory activity thatpropagates through the network to other brain regions.Alterations of oscillatory activity can be a sign of a braindisorder, and they are thought to be due to changes atthe level of tissue and local/global connectivity. Due toits ability to capture such oscillatory activity, EEG is commonly used in research and clinical fields to studythe neurophysiological bases of brain disorders, helpingdiagnosis and treatment (IV, 2014). Physiologically andtheoretically inspired computational models are able toreproduce EEG signals, offering a unique tool whichcomplements experimental approaches. The applicationof computational models reveals disease mechanisms,helps testing new clinical hypotheses, and to explorenew surgical strategies in silico . This section presentscomputational models of EEG that have been employedto study different states of consciousness - wakefulness,deep sleep, anesthesia, and disorders of consciousness- as well as diseases such as neuropsychiatric disordersand epilepsy.4.1 States of consciousnessA variety of models have been employed to investigatethe brain dynamics in different physiological brain states,such as wakefulness and deep sleep (non-rapid eye move-ment, NREM) (Hill and Tononi, 2005; Cona et al., 2014;Robinson et al., 2002; Roberts and Robinson, 2012), andpharmacological conditions, such as anesthesia (Chingand Brown, 2014; Ching et al., 2010b; Sheeba et al.,2008; Hutt and Longtin, 2010; Liley et al., 2010). Othermodeling approaches seek to elucidate the neurophys-iological mechanisms underlying the presence or theabsence of consciousness in wakefulness, NREM sleepand anesthesia, and they have crucial implications forthe study of disorders of consciousness (DoC). DoC referto a class of clinical conditions that may follow a se-vere brain injury (hypoxic/ischemic or traumatic braininjury) and include coma, vegetative state or unrespon-sive wakefulness syndrome (VS/UWS), and minimallyconscious state (MCS). Coma has been defined as astate of unresponsiveness characterized by the absenceof arousal (patients lie with their eyes closed) and, hence,of awareness. VS/UWS denotes a condition of wakeful-ness with reflex movements and without behaviouralsigns of awareness, while patients in MCS show unequiv-ocal signs of interaction with the environment.The current gold standard for clinical assessment ofconsciousness relies on the Coma Recovery Scale Re-vised (Giacino et al., 2004), which scores the abilityof patients to behaviourally respond to sensory stimulior commands. However, behavioral-based clinical diag-noses can lead to misclassification of MCS as VS/UWSbecause some patients may regain consciousness with-out recovering their ability to understand, move andcommunicate (Childs et al., 1993; Andrews et al., 1996;Schnakers et al., 2009). A great effort has been devotedto develop advanced imaging and neurophysiological10echniques for assessing covert consciousness and to im-prove diagnostic and prognostic accuracy (Edlow et al.,2017; Bodart et al., 2017; Stender et al., 2014; Brunoet al., 2011; Owen and Coleman, 2008; Stender et al.,2016). A novel neurophysiological approach to unravelthe capacity of the brain to sustain consciousness ex-ploits Transcranial Magnetic Stimulation (TMS) in com-bination with EEG (Rosanova et al., 2018; Casarottoet al., 2016). Specifically, the EEG response evoked byTMS in conscious subjects exhibits complex patterns ofactivation resulting from preserved cortical interactions.In contrast, when unconscious patients are stimulatedwith TMS, the evoked-response shows a local patternof activation, similar to the one observed in healthycontrols during NREM sleep and anesthesia.The perturbational complexity index (PCI) (Casaliet al., 2013) captures the dynamical complexity of TMS-evoked EEG potentials by means of the Lempel-Zivcompression algorithm, showing high values (low com-pressibility) for complex chains of activation typical ofthe awake state, and low values (high compressibility) forstereotypical patterns of activation typical of sleep andanesthesia. PCI has been validated on a benchmark pop-ulation of 150 conscious and unconscious controls andtested on 81 severely brain-injured patients (Casarottoet al., 2016), showing an unprecedented high sensitivity(94.7%) in discriminating conscious from unconsciousstates.A recently published modeling approach (Bensaidet al., 2019) investigates the physiological mechanismsunderlying the generation of complex or stereotypicalTMS-evoked EEG responses. The proposed brain net-work model, named COALIA, describes local dynamicsas neural masses (Wendling et al., 2002) that includepopulations of pyramidal neurons and three differenttypes of interneurons. Each neural mass describes thelocal field activity of one of 66 cortical brain regions(Desikan et al., 2006). Neural masses are connected witheach other through long-range white matter fibers asdescribed above (section 3.3). EEG signals are then sim-ulated as neural mass activity. A systematic comparisonof the complexity of simulated and real TMS-evokedEEG potentials through PCI suggested that the rhyth-mically patterned thalamocortical activity, typical ofsleep, plays a key role in disrupting the complex pat-terns of activation evoked by TMS (Bensaid et al., 2019).Indeed, this rhythmical thalamocortical activity resultsin inhibition within the cortex that prevents informationfrom propagating from one brain region to another, andthus disrupts functional integration, i.e. the ability ofthe brain to integrate information originating from dif-ferent brain regions or groups of brain regions (Tononi,1998). Functional integration is necessary, along with functional segregation, i.e. the ability of brain regions orgroups of brain regions to fulfill a certain function dis-tinct from other areas of the brain (Lord et al., 2017), togenerate complex time-varying patterns of coordinatedcortical activity that are typical of the awake brain, andthought to sustain consciousness and cognitive functions(Casali et al., 2013; Demertzi et al., 2019).4.2 Neuropsychiatric disordersDisruption of integration and segregation balance, whichis fundamental for consciousness as mentioned in sec-tion 4.1, have been linked also to several neuropsychi-atric disorders as a result of altered structural and func-tional connectivity (Bassett and Bullmore, 2009; Fornitoet al., 2015; Menon, 2011; Deco et al., 2015). Among neu-ropsychiatric disorders, as reviewed in Lord et al. (2017),Alzheimer’s disease is characterized by a decrease inlong-range functional connectivity, directly affecting in-tegration between functional modules of the brain (Stamet al., 2007; Sanz-Arigita et al., 2010). Schizophrenia hasbeen linked to a “subtle randomization” of global func-tional connectivity, such that the so-called “small-world”character of the network is disrupted (Alexander-Blochet al., 2010; Lynall et al., 2010); a small-world networkis characterized by short path lengths and strong modu-larity, network properties that are thought to promoteinformation processing in the brain (Bassett and Bull-more, 2006) (but see Hilgetag and Goulas (2016)). Lossof integration has also been observed in schizophrenia(Damaraju et al., 2014).As explained in section 3.3, whole-brain computa-tional models provide insights into how anatomical con-nections shape and constrain functional connectivity(Deco et al., 2014a, 2013a; Honey et al., 2009). UsingBNMs, Cabral and colleagues have shown that the al-terations reported in schizophrenia (Lynall et al., 2010)can be explained by a decrease in connectivity betweenbrain areas, occurring either at the local or global leveland encompassing either axonal or synaptic mechanisms,hence reinforcing the idea of schizophrenia being thebehavioural consequence of a multitude of causes dis-rupting connectivity between brain areas (Cabral et al.,2012b,a).However, these models have focused on reproducingfMRI findings and are yet to be extended to addressalterations in EEG spectral signatures in schizophre-nia, namely increased EEG gamma-band power and de-creased alpha power (Uhlhaas and Singer, 2013), which,following previous modeling insights (Cabral et al., 2014),may arise from reduced coupling between local gamma-band oscillators. Furthermore, BNMs can be employedto test how clinical interventions may help to re-establish11ealthy network properties such as balance between inte-gration and segregation or small-worldness (Deco et al.,2018, 2019a).4.3 EpilepsyModels have been employed to study pathological alter-ations of network oscillatory activity related to manydiseases, including epilepsy (Stefanescu et al., 2012;Wendling, 2005; Lytton, 2008; Holt and Netoff, 2013).Epilepsy is a complex disease which impacts 1% of theworld population and is drug resistant in approximately30% of cases. Due to its intrinsic complexity, epilepsyresearch has strongly benefited, and will do so evenmore in the future, from an in silico environment wherehypotheses about brain mechanisms of epileptic seizurescan be tested in order to guide strategies for surgical,pharmacological and electrical stimulation techniques.Focal epilepsy is a prototypical example of a diseasethat involves both local tissue and network properties.Focal epilepsy occurs when seizures originate in oneor multiple sites, so-called epileptogenic zones (EZ),before recruiting close and distal non-epileptogenic areaspertaining to the pathological network. Patients with ahistory of drug-resistant focal epilepsy are candidates forsurgery which targets epileptogenic areas and/or criticalnodes presumably involved in the epileptogenic network.Successful outcomes of these procedures critically relyon the ability of clinicians to precisely identify the EZ.A promising modeling approach aims at studying fo-cal epilepsy through a single-subject virtual brain (Proixet al., 2017; Terry et al., 2012; Hutchings et al., 2015;Bansal et al., 2018; Soltesz and Staley, 2011), bring-ing together the description of how seizures start andend (seizure onset and offset, respectively) at a locallevel (through neural mass models) (Robinson et al.,2002; Wendling et al., 2002; Lopes da Silva et al., 2003;Breakspear et al., 2006; Jirsa et al., 2014) with indi-vidual brain connectivity derived from dMRI data. Inthis personalized approach, a patient’s brain is virtu-ally reconstructed, such that systematic testing of manysurgical scenarios is possible. The individual virtualbrain approach provides clinicians with additional in-formation, helping them to identify locations which areresponsible for starting or propagating the seizure andwhose removal would therefore lead to the patient beingseizure-free while avoiding functional side effects of re-moving brain regions and connections (Olmi et al., 2019;Proix et al., 2017).
In this paper we introduced different computationalmodel types and their application to EEG, using asimple classification by spatial scale. Clearly not allmodels in the literature would necessarily belong toone category, but we believe this taxonomy can providean entry point for non-experts. The main motivationbehind this review was to identify obstacles that standin the way of applying EEG modeling in both a researchand clinical context, and to point out future directionsthat could remove these obstacles.We have pointed out several recent efforts that havebegun to more closely align models and experimentalfindings. Such integration of theory and experiment guar-antees the use of biologically relevant measures withincomputational models of EEG, a crucial element if onewishes to use EEG models together with acquired data.For example, recent microcircuit models address thegap between theory and experiment by linking averagefiring rate - a measure of population activity preferredby the modeling community - and local field potential(LFP) - a measure that is generally thought to be agood proxy of the EEG signal (Saponati et al., 2019);recent mean field models explicitly include the contri-bution of neural synchronization to the LFP (Byrneet al., 2020), thereby integrating experimental knowl-edge about how the EEG signal is generated (da Silva,2013); brain network models explore the contribution ofempirically measured connectomes to macroscopic brainactivity (Cabral et al., 2014); and applications of compu-tational models already exist that use clinical measuresto study e.g. coma (Bensaid et al., 2019), epilepsy (Olmiet al., 2019; Proix et al., 2017), and neuropsychiatricdisorders (Spiegler et al., 2016; Kunze et al., 2016). Fur-thermore, some modeling approaches focus on providinga simple model for large-scale dynamics, making resultsmore interpretable both from a theoretical and clinicalstandpoint (Abdelnour et al., 2018; Raj et al., 2020).We have reviewed computational models on threespatial scales (Figure 2 A). Each scale models qualita-tively different biological processes which can be mea-sured using distinct recording techniques (Varela et al.,2001) (Figure 1). While EEG records activity at themacro-scale, mechanisms at each scale have an impacton the EEG signal and should therefore inform its inter-pretation. Therefore, ideally, scales should be combinedto provide a complete picture of neural mechanismsunderlying EEG activity, something that started to beexplored for example in the simulation platform The Vir-tual Brain (TVB) (Sanz Leon et al., 2013; Falcon et al.,2016) or in studies showing the theoretical relationshipbetween spiking networks and mean field formulations12 able 1
Some terminology used in this paper.Functional connectivity (FC) Statistical dependencies between time series recorded from different brain regionsor simulated at different nodes. Such dependencies are taken to indicate a functionalrelatedness of the brain regions/nodes. Many measures are available, for example corre-lation between amplitude envelopes, phase locking value, imaginary coherence, etc. Seefor example Colclough et al. (2016) for an overview. Note that FC does not establish acausal relationship (Friston, 2011).Structural connectivity (SC) Also known as neuroanatomical, anatomical, or white matter connectivity. Diffusion-weighted MRI (dMRI) is able to measure the diffusion of water through brain tissue(Basser et al., 1994). As water diffuses preferably along axons rather than across theirwalls, the orientation of large fiber bundles can be inferred from dMRI via algorithmsknown as fiber tracking (Jones, 2010). Note that SC does not take into account localanatomical connections made within the gray matter, and that fiber counts or densitiesdo not allow making conclusions about the weight of that connection (Jeurissen et al.,2019). Furthermore, fiber tracking algorithms are unable to resolve ambiguities introducesby crossing fibers, and it is difficult to track long fibers.Graph A brain network model, which consists of nodes and edges (Figure 2 A, right), can beformalized as a graph (Bassett and Sporns, 2017; Sporns, 2018). This can be visualizedusing so-called adjacency matrices, which contain a weight for the edge between eachpair of nodes ((Figure 2 B)). In this sense, both FC and SC matrices are adjacencymatrices. This formalization opens up the analysis of brain networks to the tools ofgraph analysis. These tools allow for example the characterization of the graph/networkusing many different quantitative measures (Rubinov and Sporns, 2010), partitioningthe graph/network into subnetworks or modules (Bassett and Sporns, 2017; Donetti andMunoz, 2004), or classifying nodes depending on their role in the network (Hagmannet al., 2008).Random walk A random walk is a random process taking place on the graph in which a “walker” isinitiated at a node and proceeds to another node following existing edges. Edges areselected by the walker with a probability proportional to their weight. Such a simulationcan be used to approximate the dynamics of spreading activation, and enables theresearcher to approximate for example the probability that activity will spread from node i to node j given the edges that exist between them, or the time that it will take foractivity to spread from node i to node j .Laplacian The Laplace operator is ubiquitous in many physical systems and is used to describestanding waves, diffusion, heat dispersion, and many other phenomena. For a network,the Laplacian is obtained directly from the adjacency matrix (see above). An intuitiveinterpretation is that it describes the “flow” of activity along the edges.Eigenmodes Many physical systems that consist of interacting elements can vibrate at certain frequen-cies, for example the string of a violin or the vibrating sheets of Chladni (Chladni, 1802).Each system has its own set of frequencies at which it can vibrate, determined for exampleby its shape. Mathematically, these eigenmodes are obtained via eigendecomposition ofthe Laplacian (see above). (Wong and Wang, 2006; Deco et al., 2013b; Coombesand Byrne, 2019; Byrne et al., 2020). Using models inthis hierarchical manner is the only way to disentangledifferent contributions to the EEG signal without usinginvasive techniques, i.e., to distinguish neural signals(Michel and Murray, 2012; Seeber et al., 2019), volumeconduction (Hindriks et al., 2017; Lind´en et al., 2011;M¨aki-Marttunen et al., 2019; Skaar et al., 2019; Einevollet al., 2013; Tele´nczuk et al., 2017; Claude Bdard, 2009),and noise. Furthermore, brain disorders can impact brainstructure and function on any scale. Using models onmultiple scales is necessary if one wishes to understandhow pathological changes manifest in clinically measur-able EEG signals. Such an understanding would also allow to use EEG to evaluate clinical interventions thataffect the micro- or mesoscale (e.g., drugs).Models can thus play an important role as a “bridge”that connects different fields. In translational applica-tions, knowledge from basic research can be integratedinto a model and the model can be designed in sucha way that it is useful for a clinical application. Anexample for a successful “bridge” is the case of BrainComputer Interfaces. In order to realize multi-scale mod-els, researchers working on animal recordings and re-searchers focusing on non-invasive recordings in humanshave to come together with modeling experts that canincorporate findings from both fields in a model.As an outlook, EEG modeling could play an impor-tant part in future endeavors towards precision medicine,13r “personal health”. Individual brain models could beused to integrate different sources of data (EEG, fMRI,ECG, etc.) in a “virtual patient”. This could complementdata-driven approaches like connectome fingerprinting(in which individuals are identified using their individualconnectome (Finn et al., 2015; Pallar´es et al., 2018; Ab-bas et al., 2020)). The ultimate goal would be to use thisvirtual patient to tailor diagnosis and therapies aroundthe needs of the patient (Wium-Andersen et al., 2017),reducing the economical burden and patient discomfortof clinical analyses and hospitalisation. Acknowledgements
The authors would like to thank AnaHernando Ariza for the first image contained in this paper andthe creativity, time and efforts she put into transforming theimportant messages of this article in graphics. The authors areparticularly grateful to Prof. J´er´emie Lefebvre and Prof. MicahM. Murray for their contributions, guidance and mentoringduring the preparation of this article.
Funding
K.G. was funded by Schweizerischer Nationalfonds zurF¨orderung der Wissenschaftlichen Forschung Award ID:170873. J.C. was funded by the Portuguese Foundationfor Science and Technology (FCT) CEECIND/03325/2017and by projects UID/Multi/50026 (FCT), NORTE-01-0145-FEDER- 000013 and NORTE-01-0145-FEDER-000023. A.C. was supported by the Tiny Blue Dot Foun-dation. A.M. was supported by PREVIEW - Bandosalute 2018 - Tuscany Regional Government. A.R. wassupported by NIH grants R01NS092802, RF1AG062196and R56AG064873. B.F. received financial support forthis work by the Fondation Asile des aveugles (grantnumber 232933), a grantor advised by Carigest SA (num-ber 232920).
Conflict of interest
The authors declare that they have no conflict of interest.
Authors contribution
K.G. and B.F. conceived the review. J.C. and A.R. con-tributed to section 2.2. A.M. contributed to section 2.1.A.C. entirely conceived section 3. K.G. and B.F. out-lined all remaining contents and all authors contributedto the final draft and review.
References θ rhythms using biophysicallocal field potential (lfp) models. Eneuro 5(4)Chiken S, Nambu A (2016) Mechanism of deep brainstimulation: inhibition, excitation, or disruption? TheNeuroscientist 22(3):313–322Childs NL, Mercer WN, Childs HW (1993) Accuracyof diagnosis of persistent vegetative state. Neurology43(8):1465–1467, DOI 10.1212/wnl.43.8.1465Ching S, Brown EN (2014) Modeling the dynam-ical effects of anesthesia on brain circuits. Cur-rent Opinion in Neurobiology 25:116–122, DOI10.1016/j.conb.2013.12.011Ching S, Cimenser A, Purdon PL, Brown EN, KopellNJ (2010a) Thalamocortical model for a propofol-induced αα