Comparison of Causality Network Estimation in the Sensor and Source Space: Simulation and Application on EEG

The usage of methods for the estimation of the true underlying connectivity among the observed variables of a system is increasing, especially in the domain of neuroscience. Granger causality and similar concepts are employed for the estimation of the brain network from electroencephalogram (EEG) data. Also source localization techniques, such as the standardized low resolution electromagnetic tomography (sLORETA), are widely used for obtaining more reliable data in the source space. In this work, connectivity structures are estimated in the sensor and in the source space making use of the sLORETA transformation for simulated and for EEG data with episodes of spontaneous epileptiform discharges (ED). From the comparative simulation study on high-dimensional coupled stochastic and deterministic systems originating in the sensor space, we conclude that the structure of the estimated causality networks differs in the sensor space and in the source space. Moreover, different network types, such as random, small-world and scale-free, can be better discriminated on the basis of the data in the original sensor space than on the transformed data in the source space. Similarly, in EEG epochs containing epileptiform discharges, the discriminative ability of network topological indices was significantly better in the sensor compared to the source level. In conclusion, causality networks constructed at the sensor and source level, for both simulated and empirical data, exhibit significant structural differences. These observations indicate that further studies are warranted in order to clarify the exact relationship between data registered in the sensor and source space.

1 Introduction

There is an increasing interest in neuroscience in working on the true current distribution of the sources in the grey matter of the brain, termed source space, and not on the extracranial electroencephalograms (EEG) or magnetic encephalogram (MEG) recordings, termed sensor space (De Vico Fallani et al., 2010; Lehmann et al., 2014). Many strategies have been proposed in order to obtain the primary current distribution of the sources in the brain from EEG or MEG. These are basically head models that take into account the volume conduction of the brain and other properties in order to estimate the activity of the initial sources that explain best the extracranial electric potential measurements. Alternatively, it has also been proposed to reconstruct the electrostatic field within a predefined cubic grid of center points, based on the assumption of an ellipsoid and electromagnetically homogeneous head model, without using dipole modeling or other priors, an approach called 3D vector field tomography (Papadaniil and Hadjileontiadis, 2015; Papadaniil et al., 2016). The methods of source reconstruction are separated in two main classes, the over-determined and the under-determined inverse models (Michel et al., 2004). The over-determined inverse models presuppose that a small number of points in the source space is capable of explaining the extracranial measurements sufficiently. In this case, a unique solution, in terms of source location and current activity, is provided when the number of parameters to be estimated is less or equal to the number of sensor space channels. On the contrary, the under-determined inverse models consider a dense three dimensional (3D) grid of points in the brain having fixed positions and being many more than the sensor space channels, which leads to infinite solutions, as stated first in (Helmholtz, 1853). The objective for these models then is to determine a unique optimal source electric activity distribution over the grid of points.

In recent years, the concepts of Granger causality (Granger, 1969) and causality networks are of increasing interest in many branches of science such as finance (Hong et al., 2009; Billio et al., 2012), socioeconomics (Bollen et al., 2011), climatology (Runge et al., 2015; Dijkstra et al., 2019) and neuroscience (Kugiumtzis and Kimiskidis, 2015; Porta and Faes, 2016; He et al., 2019; Rossini et al., 2019). Intuitively, given two variables, X₁ and X₂, the Granger causality from variable X₁ to variable X₂ exists if the past and present of X₁ provide information that improves forecasts for the future of X₂. The concept of Granger causality is implemented in a number of measures in the time, frequency and phase domain, and extended to account for the presence of other observed variables and thus estimate only direct causal effects. In several studies, linear and nonlinear as well as bivariate and multivariate causality measures are compared, e.g., see the recent comparative studies in (Bakhshayesh et al., 2019; Siggiridou et al., 2019; Sommariva et al., 2019) and references therein, and it turns out that the most appropriate measures are the direct and nonlinear measures, but these are also the harder to be estimated reliably. The dimension reduction, intrinsically taking place in the algorithm of some Granger causality measures, allows for the estimation of direct causal effects in high-dimensional time series, such as multi-channel EEG. We consider a linear and a nonlinear such measure in this study, the restricted conditional Granger causality index (RCGCI) (Siggiridou and Kugiumtzis, 2016) and the partial mutual information from mixed embedding (PMIME) (Kugiumtzis, 2013), respectively, both found to perform well in high-dimensional settings (Koutlis and Kugiumtzis, 2016; Siggiridou et al., 2019).

The direct measures are particularly relevant for the formation of causality networks from multivariate time series because only direct causal effects are estimated. The causality networks are graphs with nodes that represent the observed variables of a system and connections between nodes that are weighted by the values of the selected causality measure. These networks are usually analyzed with simple metrics of network theory in order to estimate important characteristics of their topology (Rubinov and Sporns, 2010; Fornito et al., 2016; Geier and Lehnertz, 2017). In neuroscience, the estimation of causality effects among brain areas using EEG recordings is a widely used approach for the observation of brain reactions to certain stimuli and for the discrimination of normal and abnormal states of brain function (Schelter et al., 2006; Lehnertz et al., 2014; Fornito et al., 2016).

For the computation of the connectivity in the source space, in some cases the cortical activity may be estimated for a limited number of nodes and that makes the computations cost-efficient (De Vico Fallani et al., 2010), while in other cases, the number of nodes is prohibitive for a multivariate analysis of the estimated cortical activity (Milz et al., 2014). In the latter cases, the usual strategy is to discretize the cortex in different regions of interest (ROIs) and to consider as new nodes the centers of these regions and assign to them the average activity of all nodes of this region (Hata et al., 2016; Toppi et al., 2016). The connectivity of the estimated cortical activity is estimated using the measures of Granger causality discussed above [see also (Lei et al., 2015)], as well as linear measures (Schoffelen and Gross, 2009; De Vico Fallani et al., 2010; Lehmann et al., 2014; Milz et al., 2014) and nonlinear measures (Mulert et al., 2011; Pascual-Marqui et al., 2014) of similarity and synchronization, i.e. phase synchronization or correlation in the time and frequency domain. Recent studies found that the choice of the inverse method impacts the brain connectivity analysis, e.g., see source-space coherence analysis on magnetoencephalogram (MEG) (Hincapié et al., 2017).

To the best of our knowledge there have not been any reports for systematic comparison of the brain connectivity estimated in the sensor and source space. These two approaches are tested in simulation examples and discussed in (Brunner et al., 2016; Van de Steen et al., 2019). Furthermore, the effect of noise and source locations on the estimation of connectivity in source space is studied in (Anzolin et al., 2019). It is well accepted that EEG sensors do not capture signals exclusively from areas beneath EEG electrodes and the location and orientation of the sources influence critically the signals registered at the sensor level (Van de Steen et al., 2019). For that very reason, in the present study, we examine the overall connectivity and not distinguish one-to-one relationship between connections of individual anatomic areas in the sensor and source space, as done in the two studies in (Brunner et al., 2016; Van de Steen et al., 2019). Our study differs from the two studies also in that it assumes the simulated systems in the sensor rather than the source space, where the latter is physiologically more appropriate as sensors do not interact with each other but brain sources do. However, the data at hand are the scalp EEG (sensor space) and it is thus reasonable to assume simulation systems for the observed dynamics. Moreover, the study is on the overall structure of causality networks and not causal relationships of distinct sources as in (Anzolin et al., 2019; Van de Steen et al., 2019), and for this to develop a high dimensional system in the source space under realistic assumptions is a hard task. For the comparison of causality networks in sensor and source space it was reported strong correlation of the global functional connectivity between the two domains on real scalp EEG of eyes open and closed (Lai et al., 2018). We extend this study using effective connectivity measures and a different real EEG setting of changing connectivity structure.

In this work, the objective is beyond investigating differences in the estimated connectivity networks in the sensor and source space, which is expected due to the transform and we show it here analytically. Rather, the objective is to apply the same procedures for causality network estimation in the sensor and source space and to compare the causality structures that arise on each of the two workspaces. Though it is not expected to find similar causality structures, an important question is whether the estimated causality structures hold the same discriminative information in any of their characteristics, estimated by network metrics. For this, we performed a simulation study on systems of linear and nonlinear dynamics, having a predefined connectivity structure of one of the three types of random (RAND) (Erdös and Rényi, 1959), scale-free (SCF) (Barabási and Albert, 1999) and small-world (SW) (Watts and Strogatz, 1998) networks. Examples of the network types are presented in Figure 1. The original space of the generated multivariate time series is assumed to be the sensor space and the data are transformed by sLORETA to source space. The objective of the simulation study is the discrimination of different coupling structures (the three network types) by a network characteristic in both spaces. As established in our previous works, the two causality measures RCGCI and PMIME have a very good discrimination capability in the sensor space (the generated time series in the simulation study) and the interest is merely to investigate the deformation of the network structure in the source space. For the neuroscience data analysis, EEG data with episodes of epileptiform discharges (ED) are used in order to compare the causality structure of the brain before, during and after an ED in the sensor and source space.

FIGURE 1

FIGURE 1. Illustration of network types: (A) random, (B) small world, and (C) scale-free.

The paper is organized as follows: in Section 2 the source localization method sLORETA, the causality measures and the network indices are presented, in Section 3 the simulated and the EEG data are briefly discussed, in Section 4 the results are presented and in Section 5 conclusions are given.

2 Methods

2.1 sLORETA

A commonly used method for the 3D localization of the sources of the electromagnetic activity of the brain is the so-called low resolution electromagnetic tomography, implemented in the LORETA software (Pascual-Marqui et al., 1994). LORETA is suggested as an improvement of the minimum norm estimate (Hämäläinen and Ilmoniemi, 1994). The minimum norm estimate has a main disadvantage in that superficial sources are preferred over deep sources due to the fact that the current density vector is not weighted for the solution of the minimization problem. LORETA overcomes this problem using as weights the norm of the columns of the lead field matrix (the matrix that transforms current density of the sources to extracranial measurements), thus giving truly 3D solutions. Additionally, it provides “smooth” solutions in terms of the minimum weighted Laplacian, which are more plausible as neighboring neurons tend to have coherent waveforms. The way to achieve a meaningful 3D unique solution is by sacrificing spatial resolution, namely the brain volume is discretized in v = 2394 voxels considered as separate dipoles and thus the name low resolution electromagnetic tomography. The forward model is defined as:

$\mathbf{m}=\mathbf{G}\mathbf{d}(1)$

where $\mathbf{m}\in {\mathbb{R}}^c$ is the vector of measurements at a certain time point, $c\in \mathbb{N}$ is the number of channels, $\mathbf{d}\in {\mathbb{R}}^{3v}$ is the vector of current densities of the v voxels (for each voxel a 3D vector is considered) and $\mathbf{G}\in {\mathbb{R}}^{c\times 3v}$ is the lead field matrix. The discrete problem is:

${\text{min}}_{\mathbf{d}}\Vert \mathbf{B}\mathbf{W}\mathbf{d}\Vert \text{,\hspace{0.17em}under\hspace{0.17em}the\hspace{0.17em}constraint\hspace{0.17em}}\mathbf{m}=\mathbf{G}\mathbf{d}(2)$

where $\mathbf{W}\in {\mathbb{R}}^{3v\times 3v}$ is a diagonal matrix with W_ii = ‖G_i‖, G_i is the ith column of G and B is the discrete Laplacian operator 3v × 3v matrix [for more details on the computation of B see the Appendix of (Pascual-Marqui et al., 1994)]. The unique solution for the current densities is:

${\mathbf{d}}^{*}=\mathbf{T}\mathbf{m}(3)$

with T=(WB^TBW)⁻¹G^T [G(WB^TBW)⁻¹G^T]^†, where the † exponent denotes the Moore-Penrose pseudo-inverse.

In (Pascual-Marqui, 2002) an improvement of LORETA is presented suggesting that it achieves zero localization error due to a standardization of the current densities. This improvement is termed standardized low resolution brain electromagnetic tomography (sLORETA) and provides a better resolution with v = 6239 voxels. The followed procedure is exactly the same as in LORETA but the estimation d* is finally standardized by its variance. This is actually similar to the work conducted in (Dale et al., 2000) utilizing standardization as well but in a different way from that in sLORETA, which conversely achieves exact localization. More precisely, the covariance matrix $S_{{\mathbf{d}}^{*}}$ of the estimated current density d* is computed and the standardized current densities are:

${{\mathbf{d}}_i^{*}}^T{\left\{{[S_{{\mathbf{d}}^{*}}]}_{ii}\right\}}^{-1}{\mathbf{d}}_i^{*}(4)$

where ${\mathbf{d}}_i^{*}\in {\mathbb{R}}^3$ is the current density estimate at the ith voxel and ${[S_{{\mathbf{d}}^{*}}]}_{ii}\in {\mathbb{R}}^{3\times 3}$ is the ith diagonal block of matrix $S_{{\mathbf{d}}^{*}}$.

2.2 Causality Network Estimation

For the estimation of the causality networks we use two Granger causality measures, one linear, and one nonlinear, and both in the time domain making use of dimension reduction and being capable of estimating direct causal effects from high-dimensional time series (Koutlis and Kugiumtzis, 2016; Siggiridou et al., 2019).

2.2.1 RCGCI

The measure restricted conditional Granger causality index (RCGCI) is the well-known linear stochastic model-based measure of conditional Granger causality index (CGCI) but based on a restricted (sparse) vector autoregressive (VAR) model (Siggiridou and Kugiumtzis, 2016). The CGCI makes use of the unrestricted VAR model (U-model) of all the lagged variables up to a maximum lag (order) p so that in total there are Kp terms, where K is the number of observed variables. Applying a selection scheme using augmented VAR models, called modified backward-in-time selection, a small subset of lagged variables are selected in the sparse VAR model to constitute the U-model for RCGCI. For high-dimensional systems with sparse coupling structure, the selected subset may have cardinality much smaller than Kp. The restricted model (R-model) is derived by excluding the lagged terms of the driving variable X, so that if the sparse VAR U-model contains no lagged terms of X then RCGCI = 0. If it contains lagged terms of X then RCGCI is computed by the log ratio of the fitted error variances of the U-model and R-model, s_U² and s_R², respectively, as for the CGCI

${\text{RCGCI}}_{X\to Y}=\log \frac{s_R^2}{s_U^2}.(5)$

2.2.2 PMIME

The measure partial mutual information from mixed embedding (PMIME) is a normalized version of the partial transfer entropy (Schreiber, 2000; Papana et al., 2012) restricted to the selected mixed embedding vector, formed by selecting progressively from the set of all delayed variables the most relevant ones to the future of the response variable (Kugiumtzis, 2013). The selection relies on information measures, mutual information (MI) and conditional mutual information (CMI). In particular, in the first step of the selection scheme, the lagged variable w with the highest mutual information (MI) to the response Y (one time step ahead, y_t+1), I(y_t+1; w), is selected and constitutes the current subset of lagged variables w_t. For each subsequent step, the lagged variable w with the highest CMI to the response given the current subset, I(y_t+1; w|w_t), is selected and added to w_t. The significance of CMI (MI for the first step) is tested at each step and the selection scheme terminates when it is found statistically not significant. Similarly to RCGCI, if the derived subset does not include any lagged terms of the driving variable X then PMIME = 0. If it includes, the information of the lagged terms of X about the response is evaluated by the CMI of the lagged terms of X in w_t, ${\mathbf{w}}_t^x$, and the response given all other components in w_t, normalized by the mutual information of all the lagged terms and the response

${\text{PMIME}}_{X\to Y}=\frac{I(y_{t+1};{\mathbf{w}}_t^x|{\mathbf{w}}_t^y,{\mathbf{w}}_t^z)}{I(y_{t+1};{\mathbf{w}}_t)},(6)$

where ${\mathbf{w}}_t^y$ and ${\mathbf{w}}_t^z$ are the lagged terms of the response Y and the remaining variables Z in w_t, respectively, so that ${\mathbf{w}}_t=[{\mathbf{w}}_t^x,{\mathbf{w}}_t^y,{\mathbf{w}}_t^z]$.

Both causality measures include an internal criterion to determine the presence of the driving of one variable to another and in that case a positive value of RCGCI or PMIME quantifies the strength of the driving, otherwise the measure is zero. This leads to a weighted matrix with zero and weighted connections of driving strength, thus the binarization of the causality network is performed just by accepting the positive values as existing causal effects.

2.3 Network Metrics

We select five network metrics, which capture different characteristics of the network structure (Rubinov and Sporns, 2010). These metrics are computed on the adjacency matrices formed by the positive PMIME or RCGCI values for each ordered pair of variables and regard the causality network estimated from the multivariate time series.

The five network metrics quantify different properties of the network with binary directed or non-directed connections. The centrality property of the network is quantified by the mean, deg^m, and standard deviation (SD), deg^SD, of the degree distribution over the nodes. The characteristic path length, λ, is used as index of functional integration. The small-worldness index, SWi, quantifies the presence of small-world structure in the network, which suggests the combination of functional segregation and functional integration in a network. As an index of resilience we consider the assortativity coefficient, r, quantifying the preference of network nodes to attach to other nodes of similar degree, typically defined as the correlation coefficient between the degrees of two nodes.

2.4 Evaluation Index

In the focus of the study is the discrimination of different connectivity structures in the sensor and source space, and we consider the three network types, i.e., random (RAND), small world (SW), and scale-free (SCF) (see Figure 1). We consider a number of realizations of multivariate time series from each network type and compute the causality networks and subsequently the network metrics for each realization. For each pair formed from any two of RAND, SW, and SCF, the evaluation index area under receiver operating characteristic curve (AUROC) (Fawcett, 2006) is employed to quantify the overlapping of the two distributions of the network metric in the two network types. For the simulation study, the AUROC is computed for each pair of network types of the original causality network (given by the system equations) and the estimated causality networks in the sensor and source space.

3 Materials

In the simulation study, the synthetic systems are considered to evolve in the sensor space. One could argue to consider the source space instead, but we do not follow this line here for two reasons. First, the assumption of the sensor space as the space of the generated time series is a natural choice since the starting data, the acquired EEG measurements, are in the sensor space. Secondly, the generation of data in the source space is a hard task and can be realized under rather unrealistic simplifications, e.g., assuming the same activity in all vertices (voxels) in a ROI on the cortical mesh (Silva Pereira et al., 2017). We consider three systems of different type. The two first systems are deterministic chaotic dynamical systems while the third is a linear multivariate stochastic process. The two first systems are coupled systems, i.e., are formed from a number of subsystems of the same type coupled to each other. The first system is in discrete time and the second in continuous time.

3.1 Simulated Data

The system of coupled Hénon maps (CHM) (Kugiumtzis, 2013) is a system of coupled chaotic maps defined as

$x_{j,t}=1.4-{\left[\frac{{\sum }_{i=1}^K{C}_{ij}{x}_{i,t-1}}{{\sum }_{i=1}^K\mathrm{\Theta }(C_{ij})}+\left(1-\frac{{\sum }_{i=1}^K{C}_{ij}}{{\sum }_{i=1}^K\mathrm{\Theta }(C_{ij})}\right)x_{j,t-1}\right]}^2+0.3x_{j,t-2}(7)$

where j = 1, 2, … , K is the variable index, K denotes the number of variables and C_ij is the coupling strength (considering x_i as the driving variable and x_j as the response variable). The Θ(C_ij) is the Heaviside function, being one if C_ij > 0 and zero if C_ij = 0. For this system, we consider K = 20 variables coupled with strength C_ij = C = 0.2 for the non-vanishing terms according to each of the three network types (RAND, SW, SCF). For each network type, 50 realizations are generated and the time series have length N = 2048.

The system of coupled Mackey-Glass (CMG) (Senthilkumar et al., 2008; Kugiumtzis, 2013) is a system of coupled identical delayed differential equations defined as

${\dot{x}}_j(t)=-0.1x_j(t)+\sum _{i=1}^K\frac{C_{ij}{x}_i(t-\mathrm{\Delta })}{1+x_i{(t-\mathrm{\Delta })}^{10}},(8)$

where j = 1, 2, … , K is the variable index, and Δ is the delay parameter. For this system we consider Δ = 100 time units (for the uncoupled system each equation defines a chaotic system of correlation dimension close to 7) and C_ij = C = 0.1. The solution is run in Matlab (Natick, Massachusetts: The MathWorks Inc.) using the function “dde23” and then sampled at 4 time units. Further details for the generation of the CMG time series are given in (Kugiumtzis, 2013). The setting is as for CHM, K = 20, 50 realizations of each of the RAND, SW, and SCF network types and N = 2048.

The third systems is a linear stochastic VAR process (Basu and Michailidis, 2015) on K = 20 variables and order p = 1,

${\mathbf{x}}_t=A{\mathbf{x}}_{t-1}+{\mathbf{e}}_t,(9)$

where x_t is the state vector of length K and e_t is the white noise vector of length K following Gaussian distribution with zero mean and unity covariance matrix. The components of the square coefficient matrix A are zero or positive determined by the selected network of type RAND, SW, or SCF. Initially, the non-vanishing coefficients are set to 0.9 (the rest are set to zero) and then they are iteratively reduced until the stationarity condition is fulfilled. We get 50 multivariate time series for each of the three network types of length N = 2048.

The PMIME is used to estimate the connectivity for the two first nonlinear systems and the RCGCI for the third linear system, as being the most appropriate measures to capture the causality effects in each system. This is so because the objective here is to have an appropriate measure that estimates best the causality network from the time series in the sensor and source space.

All the computationally generated data are treated as if they were EEG signals (in the sensor space) in order to transform them with sLORETA (to the source space) and the K = 20 variables are considered as electrode positions on the scalp, more precisely: FC3, FC1, FCz, FC2, FC4, C3, C1, Cz, C2, C4, CP3, CP1, CPz, CP2, CP4, P3, P1, Pz, P2, P4. These positions are selected in order to cover a sufficient area of the brain rather than concentrating on a small region, as shown in Figure 2A.

FIGURE 2

FIGURE 2. (A) Electrode positions assumed for the sLORETA on the simulated data. (B) One epoch with an ED episode. The vertical lines represent the start and end points of the ED, respectively, and the time order is with respect to the ED start (0 time at start).

For each of the 50 realizations, the hypothesized multi-channel (K = 20) time series in the sensor space is transformed with sLORETA to v = 6239 source space signals, which consequently are separated in ROIs and averaged, thereby resulting in K = 20 source space time series. The ROIs are selected with the option of “all nearest voxels” in the sLORETA software, which separates the voxels in mutually disjoint regions of the nearest to the selected electrodes voxels.

3.2 EEG Data

For the real data analysis, EEG data from three epileptic patients at rest are used. In order to avoid the confounding variables of diverse structural substrates and varying localizations of the epileptic zone, which typically characterize focal epilepsies, we used data from patients suffering from Genetic Generalized Epilepsies. From the initial recording of 60 channels, channels containing artifacts were rejected, i.e., 16 channels for subject 1, 9 channels for subject 2, and 14 channels for subject 3. The EEG data of the artifact-free channels were band-pass filtered in (0.01,70) Hz, downsampled to 200 Hz, and re-referenced to infinity (Yao, 2001), a re-referencing scheme found to be more appropriate for connectivity analysis (Qin et al., 2010). Epochs containing epileptiform discharge (ED) were selected, 7 epochs from subject 1, 3 epochs from subject 2 and 10 epochs from subject 3. An exemplary epoch from subject 1 is shown in Figure 2B. Each epoch contains a pre-ED period of 10 s, a spontaneous ED of duration from 2.19 to 5.71 s and a post-ED period of 10 s.

The sLORETA software is used for the transformation of the EEG signals to the source space. First, the electrode positions are selected and then the transformation matrix is generated by the software. Consequently, the signals are transformed to the source space giving a large number of v = 6239 time series. In accordance with the simulation study, a number of 44 ROIs corresponding to the 44 scalp electrodes for subject 1 are selected and 44 time series for the source space are obtained, and the same was done accordingly for the artifact-free channels of the other two subjects.

After the transformation of the time series, both the sensor and source time series are split in sliding overlapping windows of 2 s with a sliding step of 0.5 s to profile the brain network characteristics and all their changes during the epochs. The causality estimation is performed with PMIME and then characteristics of the causality networks are estimated by the network metrics previously presented.

4 Results

Before presenting the results of the simulation study, it is elaborated that a linear transformation may change the connectivity of a set of variables in a significant degree. Although the employed methodology (sLORETA) is not a linear transformation, the following simplistic example is indicative of the changes a transformation in general can induce. Let us assume a linear VAR(1) process on K variables

${\mathbf{x}}_t=A{\mathbf{x}}_{t-1}+{\mathbf{e}}_t,(10)$

where the coefficient matrix A defines the connectivity matrix of the original system, x_t = [x_1,t, … , x_K,_t] is the vector of the state of K variables at indexed time t and e_t = [e_1,t, … , e_K,_t] is the vector of noise terms at time t having normal distribution e_i,_t ∼ N(0, σ²), i = 1, … , K, and covariance matrix ${\mathrm{\Sigma }}_{{\mathbf{e}}_t}={\sigma }^2{I}_K$, where I_K is the identity matrix of size K × K. Consider a linear transformation of x_t

${\mathbf{y}}_t=H{\mathbf{x}}_t,(11)$

where H is a matrix of size K × K with det(H) ≠ 0. Solving Eq. (11) with respect to x_t and substituting it in Eq. (10), the VAR (1) for y_t is obtained as:

${\mathbf{y}}_t=HA{H}^{-1}{\mathbf{y}}_{t-1}+H{\mathbf{e}}_t=B{\mathbf{y}}_{t-1}+{\mathbf{e}}_t^{\prime },(12)$

where the input noise vector ${\mathbf{e}}_t^{\prime }$ is correlated (the noise covariance matrix is not diagonal) indicating also the presence of instantaneous causality, which however is not relevant here, as we estimate only the lag causality. In analogy to the connectivity of the original VAR(1) system, the connectivity of the transformed VAR(1) system is determined by the coefficient matrix B = HAH⁻¹. Thus, a variety of different connectivity structures for the transformed system can be obtained, dependent on the form of H, with the only constraint being that the initial and the transformed system have connectivity matrices with the same eigenvalues. Also, it is noteworthy that the transformation matrix $H=Q_A^{-1}$, where the columns of Q_A are the eigenvectors of A, leads to a diagonal B and subsequently a linearly transformed VAR(1) system with no connections among its variables, regardless of its initial connectivity structure.

4.1 Simulation Study

In Figure 3, the steps of the procedure followed in the simulation study are illustrated in an example for the CMG system. Initially, a network of a predefined type (RAND, SW, or SCF) is defined, in order to form the coupling relationships in the system equations, for the generation of the 50 simulated time series of each network type in the hypothesized sensor space. Then, sLORETA is used to transform the time series to the source space. The causal effects for all pairs of observed variables are estimated by a causality measure (PMIME for CMG in Figure 3) in both the sensor and the source space, resulting in the respective causality networks. It is observed that the form of the time series is rather different before and after the sLORETA transformation. More precisely, the source space time series are less oscillatory and show spikes at various times contrary to the initial time series.

FIGURE 3

FIGURE 3. The initial scale-free coupling structure among the variables is used to generate the K = 20 time series from the CMG system hypothesized to be in the sensor space. Then, sLORETA derives the respective time series in the source space and finally the PMIME is used to estimate the causal effects and give the causality networks in both the sensor and the source levels.

Regarding the causality network estimation, for the example in Figure 3 the PMIME has a satisfactory performance on the sensor space time series capturing the original scale-free structure. The causality network of the transformed time series in the source space is qualitatively different from the original network, as expected from the previous analytical example. This is quantitatively confirmed for instance by the network metric deg^SD that takes the values 5.33 and 2.35, respectively, for the initial structure and the source space estimation. In Figure 4, aggregate results from the 50 realizations are shown for all three network structures, i.e., the initial network, the network estimated in the sensor space and the network estimated in the source space. The similarity of the deg^SD between the original network and the estimated network in the sensor space as well as their dissimilarity to the estimated network in the source space are striking for the cases of RAND and SW networks. For SCF networks, the PMIME does not estimate with the same high accuracy the initial network (see specific results in (Koutlis and Kugiumtzis, 2016)) and therefore the histograms of deg^SD for the sensor and source space networks differ significantly from the respective histogram for the initial network.

FIGURE 4

FIGURE 4. Histograms for the standard deviation of the degree, deg^SD, for 50 realizations of each network type: (A) random (RAND), (B) small-world (SW), and (C) scale-free (SCF). At each panel, there are three histograms, for the initial network, the network estimated in the sensor space, and the network estimated in the source space, as given in the legend.

The matching of the metrics of the estimated network to the true network can be seen collectively in the scatter plots for each network metric and for both RAND and SW network types in Figure 5. Results are not shown for the characteristic deg^m being constant for the true networks and the SCF network type as for all but the network metric r the metric of the true network is rather constant across the 50 realizations. The network metrics SWi, deg^SD, and λ computed on the sensor space tend to match well with the respective metrics of the true networks whereas for the source space the respective metrics are more spread and have the tendency to overestimate or underestimate the metrics of the true network (over or under the diagonal in Figures 5A–C, respectively). For the metric r in Figure 5D the values of the true networks are at the same range for RAND and SW network types so that at first look r values seem to spread about the same for both sensor and source space. A careful look would discern the alignment of the r values for the sensor space along the diagonal and indeed the Pearson correlation coefficient values are for RAND network type 0.84 for sensor space as opposed to 0.01 for source space and for SW network type 0.67 and −0.28, respectively.

FIGURE 5

FIGURE 5. Scatter plots of network metrics estimated on 50 realizations of random (RAND) and small-world (SW) network types, where in the x-axis is the metric for the true network and in the y-axis for the estimated network in the sensor and source space, as shown in the legend. The metrics are SWi in (A), deg^SD in (B), λ in (C) and r in (D).

Establishing that the connectivity changes significantly under the sLORETA transform, we investigate now whether the transform preserves the discriminative information that allows for correct classification of the different initial coupling structures. For this, we examine the discrimination of the three network types in the sensor and source space. In Figure 6, results for the discriminative ability of the network metrics are presented for the CMG system, where the estimated networks are derived by the PMIME.

FIGURE 6

FIGURE 6. Distribution of network metric values for three different network types (RAND, SW, and SCF) computed on the initial network and on the networks estimated by PMIME on the sensor space and on the source space, as shown in the legend. The PMIME is computed on simulated time series of the CMG system. The metrics are SWi in (A), deg^m in (B), deg^SD in (C), λ in (D) and r in (E).

It is observed that the values of all network metrics for the initial network and the network estimated on the sensor space range in a similar scale. This includes also the case of deg^SD, where though the initial connections are not preserved well in the sensor space with the PMIME for SCF networks (see Figure 4C), still the SCF networks can be clearly discriminated with deg^SD from RAND and SW networks. The values of the network metrics for the source space are in some cases similar with those of the initial networks, e.g., SWi on RAND network, but in most cases they range in different scales. This result indicates that the transformation to the source space changes significantly the structure of the causality network. The discriminative information of the network metrics, regarding the three network types, remains in the sensor space estimation at a larger degree than in the source space. This result is confirmed by the average AUROC values for the differences in the three pairs of the RAND, SW, and SCF network types, presented in Table 1. While the AUROC values for the networks estimated in the sensor space are high and similar to these for the initial networks, the respective AUROC values for networks estimated in the source space are much lower. However, in almost all of the discrimination tasks of Table 1 we obtain a p-value < 0.01 after ANOVA hypothesis testing for equal mean AUROC in the RAND, SW and SCF network types, with the exception of source space estimation of SWi on CMG, source space estimation of SWi on CHM, source space estimation of r on CHM and source space estimation of r on VAR.

TABLE 1

TABLE 1. Average AUROC values from the three AUROC values computed in the three binary classification tasks, RAND vs. SW, RAND vs. SCF, and SW vs. SCF, for the five network metrics (labels in the first row), the three dynamical systems (labels in the first column) and for the initial networks, the networks from the sensor space and the networks from the source space.

In Figure 7, results for the discriminative ability of the network indices are presented for the CHM system in the same way as shown earlier for the CMG system in Figure 6. Similarly to the CMG case, the values of all network indices for the initial network and the network estimated on the sensor space range in a similar scale. In contrast, the values of the network indices for the networks estimated on the source space range in different scales with those of the initial networks in most of the cases. Again, this is confirmed with the results of AUROC in Table 1.

FIGURE 7

FIGURE 7. Distribution of network metric values for three different network types (RAND, SW, and SCF) computed on the initial network and on the networks estimated by PMIME on the sensor space and on the source space, as shown in the legend. The PMIME is computed on simulated time series of the CHM system. The metrics are SWi in (A), deg^m in (B), deg^SD in (C), λ in (D) and r in (E).

In Figure 8, the respective results for the VAR system are presented. Here, the causality networks are estimated with the linear causality measure RCGCI, which is more appropriate (and simple) for the linear VAR process. Regarding the range of network index values, in many cases the initial network, the sensor space network and the source space network share the same range of values, e.g., for the random network type and all measures or the scale-free network and SWi and λ. Also, the discriminative information of the network indices is preserved in the sensor space in all cases, but in the source space only for SWi and λ, as it is shown also in Table 1.

FIGURE 8

FIGURE 8. Distribution of network metric values for three different network types (RAND, SW, and SCF) computed on the initial network and on the networks estimated by PMIME on the sensor space and on the source space, as shown in the legend. The PMIME is computed on simulated time series of the VAR system. The metrics are SWi in (A), deg^m in (B), deg^SD in (C), λ in (D) and r in (E).

4.2 EEG Data Analysis

For the EEG data analysis, the epochs from three subjects that contain epileptiform discharges (ED) are considered. The EDs of the 7 epochs of Subject 1 (S1) have duration from 2.19 to 5.71 s and the EDs of the 9 epochs of Subject 3 (S3) from 1.94 to 5.63 s, whereas the 10th epoch is exceptionally long at 29.05 s. The 3 epochs of Subject 2 (S2) have all very short EDs of duration 1.40, 1.80, and 2.09 s. Each epoch consists of a 10 s pre-ED period, an ED period and a 10 s post-ED period. The epoch is split in sliding overlapping windows of 2 s duration and 0.5 s step on which the causality estimation is performed with PMIME. As shown in Figure 9, for each window a causality network is estimated from the initial time series and from the time series obtained after the sLORETA transformation.

FIGURE 9

FIGURE 9. Upper panel: One epoch of S1 comprising a 10 s pre-ED period, an ED lasting 5.71 s and a 10 s post-ED period. Three time windows are highlighted (one in the pre-ED period, one during the ED and one in the post-ED period) and the respective causality networks are presented right under the arrows. Lower panel: the same but for the sLORETA transformed time series.

Characteristics of the causality networks are estimated by network metrics at each time window and the profiles of the metrics are shown in Figure 10 for S1, Figure 11 for S2, and Figure 12 for S3. It is observed that the network indices deg^m and λ in both the sensor and the source space have similar estimated values during the epochs, with subtle differences at various times, e.g., for S1 the coefficient of determination R² is 0.52 and 0.28, respectively. In contrast, SWi shows significant differences in the source space at many time points, e.g., for S1 R² = 0.02, which is in agreement with the simulation study.

FIGURE 10

FIGURE 10. The profile of five network indices on the sensor (green line with “o” marker) and on the source space (red line with “+” marker) for 7 epochs of S1 given in the row panels: (A)SWi, (B)deg^m, (C)deg^SD, (D)λ, (E)r.

FIGURE 11

FIGURE 11. The profile of five network indices on the sensor (green line with “o” marker) and on the source space (red line with “+” marker) for 3 epochs of S2 given in the row panels: (A)SWi, (B)deg^m, (C)deg^SD, (D)λ, (E)r.

FIGURE 12

FIGURE 12. The profile of five network indices on the sensor (green line with “o” marker) and on the source space (red line with “+” marker) for 10 epochs of S3 given in the row panels: (A)SWi, (B)deg^m, (C)deg^SD, (D)λ, (E)r.

The interest is mainly in the pre-ED and ED periods where the brain state changes abruptly at the ED onset, whereas from ED to post-ED the change is gradual as the brain state recovers from the ED. To quantify the differences in the three states, before the ED (pre-ED), within the ED (ED) and after the ED (post-ED), we computed the AUROC for pre-ED vs. ED, pre-ED vs. post-ED and ED vs. post-ED in both spaces (sensor and source) and all network indices. The results for S1 and S3 are shown in Table 2, while for S2 statistical comparison is not possible as the EDs in the three episodes are very short.

TABLE 2

TABLE 2. Average AUROC values over the 7 epochs of S1 and the 10 epochs of S3 for three discrimination tasks (pre-ED vs. ED, pre-ED vs. post-ED, ED vs. post-ED), both spaces (sensor, source) and five network metrics (labels in the first row).

The network index deg^m discriminates well the pre-ED from ED period in all episodes of the two subjects with average AUROC value 0.88 for the sensor space for both S1 and S3 and 0.80 for S1 and 0.84 for S3 for the source space. In the post-ED period, deg^m recovers its initial level gradually.

For the network index λ the respective AUROC values are 0.85 and 0.84 for S1 and S3 for the sensor space while for the source space the corresponding AUROC are 0.75 and 0.78, respectively. The other two indices, SW and r, have low discrimination ability, but still the AUROC values are higher for the sensor space than for the source space for both subjects S1 and S3. These findings are also confirmed by the independent sample t-test for equal means of pre-ED and ED, as shown in Table 3, where for more episodes at the sensor space than at the source space the discrimination is established for the stringent significance level of 0.01.

TABLE 3

TABLE 3. Results of the independent sample t-test for equal means of pre-ED and ED for subjects S1 and S3 and the five network indices as given in the first row. The numbers correspond to the ED episodes of S1 and S3 for which the p-value of the test is < 0.01 (for each subject first row is for the sensor space and second row for the source space).

The results show that the network characteristics do change after the transformation to the source space at a varying degree, so that the discrimination of the pre-ED vs. ED period is less significant. This result is not in agreement to a recent report of qualitatively similar results of connectivity structure in the source activity from reconstructed scalp EEG data and the connectivity from corresponding electrocorticographic sources in primates (Macaca mulatta) (Papadopoulou et al., 2019).

5 Discussion

In this work, the level of preservation of the main network structure when estimated on the original data, also called the sensor space, and data transformed to the so-called source space has been investigated. The transformation is performed using the software sLORETA and the causality effects among the system variables are estimated with a linear method (RCGCI) and a nonlinear information based method (PMIME). After the causality network estimation, certain network characteristics are estimated by network indices in order to compare the respective network topologies in the two spaces (sensor and source).

We performed a simulation study first in order to compare the results to the predefined ground truth, defined in the sensor space. Admittedly, the true brain system is in source space and physiologically it would be more appropriate to define the simulation systems in the source space as done in (Brunner et al., 2016; Anzolin et al., 2019; Van de Steen et al., 2019). In contrast to these studies, where causal relationships of distinct sources are investigated, our study is on the overall structure of causality networks, requiring large scale systems that could only be realized under simplifications in the source space, e.g,. in (Silva Pereira et al., 2017) this was attempted assuming the same activity in all vertices (voxels) in a ROI on the cortical mesh. On the other hand, it is common to assume the simulation systems at the domain of the acquired data (here the scalp EEG). Thus, the simulation systems refer to the observed dynamics, and the focus of the study is on the change of the overall network structure of the high-dimensional coupled system under the inverse transform defined by sLORETA.

Three different network types were considered as initial coupling structures for three high-dimensional coupled dynamical systems: K = 20 coupled Mackey-Glass systems (CMG), K = 20 Hénon coupled maps (CHM) and the linear vector autoregressive process of order one on K = 20 variables (VAR). The objective was first to designate the differences in the topology of the causality networks in the two spaces, and then to assess whether the information the network metrics hold regarding the discrimination of the three network types was preserved after the transformation to the source space. Then we proceeded to EEG data analysis with illustrative cases from recordings of three epileptic patients that contain epileptiform discharges (EDs). The objective was again to compare the differences in the topology of the estimated causality networks in the sensor and in the source space. In addition, we tried to clarify the degree to which the topology of the networks held discriminative information regarding the pre-ED period vs. the ED and in which space the discrimination was more clear.

In the simulation study, first we showed with an analytical example that a linear transformation of the system variables can change considerably the coupling structure of the dynamical system, but it can preserve some characteristics e.g. the eigenvalues of its adjacency matrix (similar derivation was obtained in Van de Steen et al., 2019). Then, the simulation results showed that the estimated causality networks have, as expected, considerable differences in terms of the topology characteristics. Finally, the ability of certain characteristics to discriminate the initial coupling structures is reduced at a varying degree after the transformation to the source space. We note here that causality estimation inaccuracy (false positives or negatives) exists also in the sensor space (the domain of true dynamics) due to many sources, such as the time series length, the size in conjunction with the density of the true coupling network, and the type of inherent dynamics and causal relationships (the latter may encounter the common source problem related to volume conduction for EEG). However, the utilized measures RCGCI and PMIME were found to perform well in high-dimensional coupled systems and therefore they were used here to allow for a better assessment of the differences in the estimated networks in sensor and source space.

In the EEG data analysis, we showed that the characteristics of the causality network topology were altered at varying degree. Here, the small-worldness index was found to change most under the transform to source space and had a very different profile at overlapping windows across the epochs that contain ED. The average degree of the causality networks was the most discriminative characteristic regarding the pre-ED vs. ED task. Also the characteristic path length showed good discrimination ability for the same task. For both indices the ability to discriminate the two periods was better in the sensor space than in the source space. The clinical EEG data are essentially a case study providing complementary evidence that clearly corroborates the basic conclusion of the simulation analysis. Further studies, on the basis of our preliminary findings, are warranted including a larger number of subjects.

We have tactically left out of the design of the connectivity analysis the issue of field spread and volume conduction. These factors may affect the results of connectivity analysis in the sensor space and subsequently the source space. It is not known as to what extent the volume conduction affects through the inverse transform the connectivity analysis in the source space (Schoffelen and Gross, 2009), but there is reported evidence from simulation studies that it does (Anzolin et al., 2019). Functional connectivity measures suggested in the literature to account for volume conduction, such as the imaginary coherence (Nolte et al., 2004), have been considered when comparing connectivity in sensor and source space (Schoffelen and Gross, 2009; Lai et al., 2018; Demuru et al., 2020). However, these measures are bivariate (allow for indirect connections) and thus not suitable to estimate the network structure, which is the main objective of this work. The utilized causality measures here are multivariate and when estimating a causal effect they account for the presence of all other observed significant sources, including to some extent the common source due to volume conduction. Certainly, they do not explicitly account for zero-lag effects brought about by volume conduction. One of the two measures used in the study has been recently modified to account for zero-lag effects, called PMIME0 (Koutlis et al., 2019). We focused on measures of lag-causality here and leave extension to other causality measures including zero-lag effects to future work.

We have also refrained from discussing the different solutions for the inverse transform and chose the sLORETA as one of the most popular transforms. Certainly, the estimated connectivity network in source space depends on the utilized inverse transform, and this is also an open research topic gaining so far little attention, e.g., two tested inverse transforms were found to give distinctly different connectivity network characteristics in source space (Lai et al., 2018).

In view of the increasing usage of source analysis in diverse areas of neuroscience, the above data suggest that further studies in order to clarify the relationship between sensor- and source-derived data are warranted.

Data Availability Statement

The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

Ethics Statement

The studies involving human participants were reviewed and approved by Aristotle University of Thessaloniki, Institutional Review Board. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

CK developed the computational procedure, made the computations and evaluation of results and wrote the manuscript. VK gave the EEG data and contributed in the discussion of the results and the editing of the manuscript. DK designed the methods and simulated systems, had the overview of the simulation study and edited the manuscript.

Funding

Part of the research was supported by a grant in the frame of the Greek–Israel Agreement on Bilateral Industrial R and D Cooperation from General Secretariat for Research and Technology of the Ministry of Education, Research, and Religious Affairs in Greece.

Conflict of Interest

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.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Anzolin, A., Presti, P., Van de Steen, F., Astolfi, L., Haufe, S., and Marinazzo, D. (2019). Quantifying the Effect of Demixing Approaches on Directed Connectivity Estimated between Reconstructed EEG Sources. Brain Topogr. 32 (4), 655–674. doi:10.1007/s10548-019-00705-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Bakhshayesh, H., Fitzgibbon, S. P., Janani, A. S., Grummett, T. S., and Pope, K. J. (2019). Detecting Connectivity in EEG: A Comparative Study of Data-Driven Effective Connectivity Measures. Comput. Biol. Med. 111, 103329. doi:10.1016/j.compbiomed.2019.103329

PubMed Abstract | CrossRef Full Text | Google Scholar

Barabási, A.-L., and Albert, R. (1999). Emergence of Scaling in Random Networks. Science 286 (5439), 509–512. doi:10.1126/science.286.5439.509

PubMed Abstract | CrossRef Full Text | Google Scholar

Basu, S., and Michailidis, G. (2015). Regularized Estimation in Sparse High-Dimensional Time Series Models. Ann. Stat. 43 (4), 1535–1567. doi:10.1214/15-aos1315

CrossRef Full Text | Google Scholar

Billio, M., Getmansky, M., Lo, A. W., and Pelizzon, L. (2012). Econometric Measures of Connectedness and Systemic Risk in the Finance and Insurance Sectors. J. Financ. Econ. 104 (3), 535–559. doi:10.1016/j.jfineco.2011.12.010

CrossRef Full Text | Google Scholar

Bollen, J., Mao, H., and Zeng, X. (2011). Twitter Mood Predicts the Stock Market. J. Comput. Sci. 2 (1), 1–8. doi:10.1016/j.jocs.2010.12.007

CrossRef Full Text | Google Scholar

Brunner, C., Billinger, M., Seeber, M., Mullen, T. R., and Makeig, S. (2016). Volume Conduction Influences Scalp-Based Connectivity Estimates. Front. Comput. Neurosci. 10, 121. doi:10.3389/fncom.2016.00121

PubMed Abstract | CrossRef Full Text | Google Scholar

Dale, A. M., Liu, A. K., Fischl, B. R., Buckner, R. L., Belliveau, J. W., Lewine, J. D., et al. (2000). Dynamic Statistical Parametric Mapping. Neuron 26, 55–67. doi:10.1016/s0896-6273(00)81138-1

PubMed Abstract | CrossRef Full Text | Google Scholar

De Vico Fallani, F., Maglione, A., Babiloni, F., Mattia, D., Astolfi, L., Vecchiato, G., et al. (2010). Cortical Network Analysis in Patients Affected by Schizophrenia. Brain Topogr. 23 (2), 214–220. doi:10.1007/s10548-010-0133-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Demuru, M., La Cava, S. M., Pani, S. M., and Fraschini, M. (2020). A Comparison between Power Spectral Density and Network Metrics: An EEG Study. Biomed. Signal Process. Control. 57, 101760. doi:10.1016/j.bspc.2019.101760

CrossRef Full Text | Google Scholar

Dijkstra, H. A., Hernández-García, E., Masoller, C., and Barreiro, M. (2019). Networks in Climate. Cambridge: Cambridge University Press.

Erdös, P., and Rényi, A. (1959). On Random Graphs I. Publicationes Math. Debrecen. 6, 290–297.

Google Scholar

Fawcett, T. (2006). An Introduction to ROC Analysis. Pattern Recognit. Lett. 27 (8), 861–874. cited By (since 1996)2503. doi:10.1016/j.patrec.2005.10.010

CrossRef Full Text | Google Scholar

Fornito, A., Zalesky, A., and Bullmore, E. T. (2016). Fundamentals of Brain Network Analysis. 1 edition. Academic Press, Elsevier.

Geier, C., and Lehnertz, K. (2017). Which Brain Regions Are Important for Seizure Dynamics in Epileptic Networks? Influence of Link Identification and EEG Recording Montage on Node Centralities. Int. J. Neural Syst. 27 (1), 1650033. doi:10.1142/s0129065716500337

PubMed Abstract | CrossRef Full Text | Google Scholar

Granger, C. W. J. (1969). Investigating Causal Relations by Econometric Models and Cross-Spectral Methods. Econometrica 37 (3), 424–438. doi:10.2307/1912791

CrossRef Full Text | Google Scholar

Hämäläinen, M. S., and Ilmoniemi, R. J. (1994). Interpreting Magnetic fields of the Brain: Minimum Norm Estimates. Med. Biol. Eng. Comput. 32 (1), 35–42. doi:10.1007/bf02512476

PubMed Abstract | CrossRef Full Text | Google Scholar

Hata, M., Kazui, H., Tanaka, T., Ishii, R., Canuet, L., Pascual-Marqui, R. D., et al. (2016). Functional Connectivity Assessed by Resting State EEG Correlates with Cognitive Decline of Alzheimer's Disease - an eLORETA Study. Clin. Neurophysiol. 127 (2), 1269–1278. doi:10.1016/j.clinph.2015.10.030

PubMed Abstract | CrossRef Full Text | Google Scholar

He, B., Astolfi, L., Valdes-Sosa, P. A., Marinazzo, D., Palva, S. O., Benar, C.-G., et al. (2019). Electrophysiological Brain Connectivity: Theory and Implementation. IEEE Trans. Biomed. Eng. 66 (7), 2115–2137. doi:10.1109/tbme.2019.2913928

PubMed Abstract | CrossRef Full Text | Google Scholar

Helmholtz, H. (1853). Ueber einige Gesetze der Vertheilung elektrischer Ströme in körperlichen Leitern mit Anwendung auf die thierisch-elektrischen Versuche. Ann. Phys. Chem. 165 (6), 211–233. doi:10.1002/andp.18531650603

CrossRef Full Text | Google Scholar

Hincapié, A.-S., Kujala, J., Mattout, J., Pascarella, A., Daligault, S., Delpuech, C., et al. (2017). The Impact of MEG Source Reconstruction Method on Source-Space Connectivity Estimation: A Comparison Between Minimum-Norm Solution and Beamforming. Neuroimage 156, 29–42. doi:10.1016/j.neuroimage.2017.04.038

PubMed Abstract | CrossRef Full Text | Google Scholar

Hong, Y., Liu, Y., and Wang, S. (2009). Granger Causality in Risk and Detection of Extreme Risk Spillover Between Financial Markets. J. Econom. 150, 271–287. doi:10.1016/j.jeconom.2008.12.013

CrossRef Full Text | Google Scholar

Koutlis, C., and Kugiumtzis, D. (2016). Discrimination of Coupling Structures Using Causality Networks from Multivariate Time Series. Chaos 26 (9), 093120. doi:10.1063/1.4963175

PubMed Abstract | CrossRef Full Text | Google Scholar

Koutlis, C., Kimiskidis, V. K., and Kugiumtzis, D. (2019). Identification of Hidden Sources by Estimating Instantaneous Causality in High-Dimensional Biomedical Time Series. Int. J. Neural Syst. 29 (4), 1850051. doi:10.1142/S012906571850051X

PubMed Abstract | CrossRef Full Text | Google Scholar

Kugiumtzis, D., and Kimiskidis, V. K. (2015). Direct Causal Networks for the Study of Transcranial Magnetic Stimulation Effects on Focal Epileptiform Discharges. Int. J. Neur. Syst. 25, 1550006. doi:10.1142/s0129065715500069

PubMed Abstract | CrossRef Full Text | Google Scholar

Kugiumtzis, D. (2013). Direct-coupling Information Measure from Nonuniform Embedding. Phys. Rev. E Stat. Nonlin Soft Matter Phys. 87, 062918. doi:10.1103/PhysRevE.87.062918

PubMed Abstract | CrossRef Full Text | Google Scholar

Lai, M., Demuru, M., Hillebrand, A., and Fraschini, M. (2018). A Comparison between Scalp- and Source-Reconstructed EEG Networks. Sci. Rep. 8, 12269. doi:10.1038/s41598-018-30869-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Lehmann, D., Faber, P. L., Pascual-Marqui, R. D., Milz, P., Herrmann, W. M., Koukkou, M., et al. (2014). Functionally Aberrant Electrophysiological Cortical Connectivities in First Episode Medication-Naive Schizophrenics from Three Psychiatry Centers. Front. Hum. Neurosci. 8, 635. doi:10.3389/fnhum.2014.00635

PubMed Abstract | CrossRef Full Text | Google Scholar

Lehnertz, K., Ansmann, G., Bialonski, S., Dickten, H., Geier, C., and Porz, S. (2014). Evolving Networks in the Human Epileptic Brain. Physica D 267, 7–15. doi:10.1016/j.physd.2013.06.009

CrossRef Full Text | Google Scholar

Lei, X., Wu, T., and Valdes-Sosa, P. A. (2015). Incorporating Priors for EEG Source Imaging and Connectivity Analysis. Front. Neurosci. 9, 284. doi:10.3389/fnins.2015.00284

PubMed Abstract | CrossRef Full Text | Google Scholar

Michel, C. M., Murray, M. M., Lantz, G., Gonzalez, S., Spinelli, L., and Grave De Peralta, R. (2004). EEG Source Imaging. Clin. Neurophysiol. 115 (10), 2195–2222. doi:10.1016/j.clinph.2004.06.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Milz, P., Faber, P. L., Lehmann, D., Kochi, K., and Pascual-Marqui, R. D. (2014). sLORETA Intracortical Lagged Coherence During Breath Counting in Meditation-Naïve Participants. Front. Hum. Neurosci. 8, 303. doi:10.3389/fnhum.2014.00303

PubMed Abstract | CrossRef Full Text | Google Scholar

Mulert, C., Kirsch, V., Pascual-Marqui, R., Mccarley, R. W., and Spencer, K. M. (2011). Long-range Synchrony of Gamma Oscillations and Auditory Hallucination Symptoms in Schizophrenia. Int. J. Psychophysiol. 79 (1), 55–63. doi:10.1016/j.ijpsycho.2010.08.004

CrossRef Full Text | Google Scholar

Nolte, G., Bai, O., Wheaton, L., Mari, Z., Vorbach, S., and Hallett, M. (2004). Identifying True Brain Interaction from EEG Data Using the Imaginary Part of Coherency. Clin. Neurophysiol. 115 (10), 2292–2307. doi:10.1016/j.clinph.2004.04.029

PubMed Abstract | CrossRef Full Text | Google Scholar

Papadaniil, C. D., and Hadjileontiadis, L. J. (2015). Tomographic Reconstruction of 3-D Irrotational Vector fields via a Discretized ray Transform. J. Math. Imaging Vis. 52 (2), 285–302. doi:10.1007/s10851-015-0559-y

CrossRef Full Text | Google Scholar

Papadaniil, C. D., Kosmidou, V. E., Tsolaki, A., Tsolaki, M., Kompatsiaris, I., and Hadjileontiadis, L. J. (2016). Cognitive MMN and P300 in Mild Cognitive Impairment and Alzheimer's Disease: A High Density EEG-3D Vector Field Tomography Approach. Brain Res. 1648 (Pt A), 425–433. doi:10.1016/j.brainres.2016.07.043

PubMed Abstract | CrossRef Full Text | Google Scholar

Papadopoulou, M., Friston, K., and Marinazzo, D. (2019). Estimating Directed Connectivity from Cortical Recordings and Reconstructed Sources. Brain Topogr. 32 (4), 741–752. doi:10.1007/s10548-015-0450-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Papana, A., Kugiumtzis, D., and Larsson, P. G. (2012). Detection of Direct Causal Effects and Application in the Analysis of Electroencephalograms from Patients with Epilepsy. Int. J. Bifurcation Chaos 22 (9), 1250222. doi:10.1142/s0218127412502227

CrossRef Full Text | Google Scholar

Pascual-Marqui, R. D., Michel, C. M., and Lehmann, D. (1994). Low Resolution Electromagnetic Tomography: A New Method for Localizing Electrical Activity in the Brain. Int. J. Psychophysiol. 18 (1), 49–65. doi:10.1016/0167-8760(84)90014-x

CrossRef Full Text | Google Scholar

Pascual-Marqui, R. D., Biscay, R. J., Bosch-Bayard, J., Lehmann, D., Kochi, K., Kinoshita, T., et al. (2014). Assessing Direct Paths of Intracortical Causal Information Flow of Oscillatory Activity with the Isolated Effective Coherence (Icoh). Front. Hum. Neurosci. 8, 448. doi:10.3389/fnhum.2014.00448

PubMed Abstract | CrossRef Full Text | Google Scholar

Pascual-Marqui, R. D. (2002). Standardized Low-Resolution Brain Electromagnetic Tomography (sLORETA): Technical Details. Methods Find Exp. Clin. Pharmacol. 24 Suppl D, 5–12.

PubMed Abstract | Google Scholar

Porta, A., and Faes, L. (2016). Wiener-Granger Causality in Network Physiology with Applications to Cardiovascular Control and Neuroscience. Proc. IEEE 104 (2), 282–309. doi:10.1109/jproc.2015.2476824

CrossRef Full Text | Google Scholar

Qin, Y., Xu, P., and Yao, D. (2010). A Comparative Study of Different References for EEG Default Mode Network: The Use of the Infinity Reference. Clin. Neurophysiol. 121 (12), 1981–1991. doi:10.1016/j.clinph.2010.03.056

PubMed Abstract | CrossRef Full Text | Google Scholar

Rossini, P. M., Di Iorio, R., Bentivoglio, M., Bertini, G., Ferreri, F., Gerloff, C., et al. (2019). Methods for Analysis of Brain Connectivity: An IFCN-Sponsored Review. Clin. Neurophysiol. 130 (10), 1833–1858. doi:10.1016/j.clinph.2019.06.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Rubinov, M., and Sporns, O. (2010). Complex Network Measures of Brain Connectivity: Uses and Interpretations. Neuroimage 52 (3), 1059–1069. doi:10.1016/j.neuroimage.2009.10.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Runge, J., Petoukhov, V., Donges, J. F., Hlinka, J., Jajcay, N., Vejmelka, M., et al. (2015). Identifying Causal Gateways and Mediators in Complex Spatio-Temporal Systems. Nat. Commun. 6, 8502. doi:10.1038/ncomms9502

PubMed Abstract | CrossRef Full Text | Google Scholar

Schelter, B., Winterhalder, M., Maiwald, T., Brandt, A., Schad, A., Schulze-Bonhage, A., et al. (2006). Testing Statistical Significance of Multivariate Time Series Analysis Techniques for Epileptic Seizure Prediction. Chaos 16 (1), 013108. doi:10.1063/1.2137623

PubMed Abstract | CrossRef Full Text | Google Scholar

Schoffelen, J.-M., and Gross, J. (2009). Source Connectivity Analysis with MEG and EEG. Hum. Brain Mapp. 30 (6), 1857–1865. doi:10.1002/hbm.20745

PubMed Abstract | CrossRef Full Text | Google Scholar

Schreiber, T. (2000). Measuring Information Transfer. Phys. Rev. Lett. 85 (2), 461–464. doi:10.1103/physrevlett.85.461

PubMed Abstract | CrossRef Full Text | Google Scholar

Senthilkumar, D. V., Lakshmanan, M., and Kurths, J. (2008). Transition from Phase to Generalized Synchronization in Time-Delay Systems. Chaos 18 (2), 023118. doi:10.1063/1.2911541

PubMed Abstract | CrossRef Full Text | Google Scholar

Siggiridou, E., and Kugiumtzis, D. (2016). Granger Causality in Multivariate Time Series Using a Time-Ordered Restricted Vector Autoregressive Model. IEEE Trans. Signal. Process. 64 (7), 1759–1773. doi:10.1109/tsp.2015.2500893

CrossRef Full Text | Google Scholar

Siggiridou, E., Koutlis, C., Tsimpiris, A., and Kugiumtzis, D. (2019). Evaluation of Granger Causality Measures for Constructing Networks from Multivariate Time Series. Entropy 21 (11), 1080. doi:10.3390/e21111080

CrossRef Full Text | Google Scholar

Silva Pereira, S., Hindriks, R., Mühlberg, S., Maris, E., van Ede, F., Griffa, A., et al. (2017). Effect of Field Spread on Resting-State Magneto Encephalography Functional Network Analysis: A Computational Modeling Study. Brain Connect. 7 (9), 541–557. doi:10.1089/brain.2017.0525

PubMed Abstract | CrossRef Full Text | Google Scholar

Sommariva, S., Sorrentino, A., Piana, M., Pizzella, V., and Marzetti, L. (2019). A Comparative Study of the Robustness of Frequency-Domain Connectivity Measures to Finite Data Length. Brain Topogr. 32 (4), 675–695. doi:10.1007/s10548-017-0609-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Toppi, J., Astolfi, L., Poudel, G. R., Innes, C. R. H., Babiloni, F., and Jones, R. D. (2016). Time-varying Effective Connectivity of the Cortical Neuroelectric Activity Associated with Behavioural Microsleeps. NeuroImage 124, 421–432. doi:10.1016/j.neuroimage.2015.08.059

PubMed Abstract | CrossRef Full Text | Google Scholar

Van de Steen, F., Faes, L., Karahan, E., Songsiri, J., Valdes-Sosa, P. A., and Marinazzo, D. (2019). Critical Comments on EEG Sensor Space Dynamical Connectivity Analysis. Brain Topogr. 32 (4), 643–654. doi:10.1007/s10548-016-0538-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Watts, D. J., and Strogatz, S. H. (1998). Collective Dynamics of ‘small-World' Networks. Nature 393 (6684), 440–442. doi:10.1038/30918

PubMed Abstract | CrossRef Full Text | Google Scholar

Yao, D. (2001). A Method to Standardize a Reference of Scalp EEG Recordings to a point at Infinity. Physiol. Meas. 22 (4), 693–711. doi:10.1088/0967-3334/22/4/305

PubMed Abstract | CrossRef Full Text | Google Scholar