Mathematical Relations Between Measures of Brain Connectivity Estimated From Electrophysiological Recordings for Gaussian Distributed Data

A large variety of methods exist to estimate brain coupling in the frequency domain from electrophysiological data measured, e.g., by EEG and MEG. Those data are to reasonable approximation, though certainly not perfectly, Gaussian distributed. This work is based on the well-known fact that for Gaussian distributed data, the cross-spectrum completely determines all statistical properties. In particular, for an infinite number of data, all normalized coupling measures at a given frequency are a function of complex coherency. However, it is largely unknown what the functional relations are. We here present those functional relations for six different measures: the weighted phase lag index, the phase lag index, the absolute value and imaginary part of the phase locking value (PLV), power envelope correlation, and power envelope correlation with correction for artifacts of volume conduction. With the exception of PLV, the final results are simple closed form formulas. In an excursion we also discuss differences between short time Fourier transformation and Hilbert transformation for estimations in the frequency domain. We tested in simulations of linear and non-linear dynamical systems and for empirical resting state EEG on sensor level to what extent a model, namely the respective function of coherency, can explain the observed couplings. For empirical data we found that for measures of phase-phase coupling deviations from the model are in general minor, while power envelope correlations systematically deviate from the model for all frequencies. For power envelope correlation with correction for artifacts of volume conduction the model cannot explain the observed couplings at all. We also analyzed power envelope correlation as a function of time and frequency in an event related experiment using a stroop reaction task and found significant event related deviations mostly in the alpha range.

1. Introduction

Electrophysiological recordings like electroencephalography (EEG) and magnetoencephalography (MEG) have a high temporal resolution, but are also non-invasive measurements with a low spatial resolution. The high temporal resolution allows to study brain oscillations, which are a ubiquitous phenomenon in many different frequency bands ranging from slow oscillations (around 1 Hz) to the high gamma rhythm (up to around 150 Hz). It is argued by many researchers that the functional role of these oscillations is a mechanism of communication between different brain areas (Engel et al., 2001, 2013; Fries, 2005, 2015). However, it is largely unclear what features of these oscillations are relevant for which specific communication within the brain.

Oscillations at a given time point, or rather segment of time, can be characterized by the frequency, the amplitude and the phase. In principle, each of these features may serve as an independent constituent of the mechanism of the communication. It is, e.g., conceivable, that phases at two neuronal sites are strongly coupled while the amplitudes are completely independent of each other and vice versa. To study the mechanisms, measures of functional dependence have been developed which mainly focus on three kinds of coupling: phase-phase-coupling, phase-amplitude coupling and amplitude-amplitude coupling (Engel et al., 2013).

The question to be addressed here is whether the corresponding measures really describe different phenomena as it is also, at least mathematically, conceivable that, e.g., phase-phase coupling determines amplitude-amplitude coupling even if the actual values are not identical, which is indeed the case for Gaussian distributed data as will be shown below. In such a case the latter would be a function of the former; the estimation of the latter would not add information on the brain dynamics and our measures would be essentially redundant. Such a redundancy occurs if the data are Gaussian distributed. In that case linear statistics, i.e., means and cross-correlation matrices or means and cross-spectra in the Fourier domain, completely determine all statistical properties. Furthermore, all coupling measures considered in this paper are normalized and independent of global (i.e., time independent) scale transformations of the data, and then all measures must be functions of complex coherency, which is the normalized version of a cross-spectrum (Nunez et al., 1997). Linear coupling measures are typically statistically robust and can be calculated with relatively low computational cost. The question then is whether the calculation of non-linear coupling is worth the effort.

In general, data are Gaussian distributed if the underlying dynamical system is linear and stationary. While, EEG and MEG data are surely not perfectly Gaussian distributed, assuming the data to be Gaussian distributed can still be a reasonable approximation. The validity of such an approximation is implicitly or explicitly assumed when estimating brain connectivity from fitting a linear dynamical model to the data as is done frequently for directed measures of connectivity like for Granger Causality (Bressler and Seth, 2011), partial directed coherence (Baccala and Sameshima, 2001), or the directed transfer function (Kaminski, 1991).

EEG and MEG have a low spatial resolution, and as consequence estimates of neuronal activities are in general mixtures of the true sources. Non-vanishing functional dependencies between such signals can be a result of such mixtures even if the underlying sources themselves are uncoupled (Nunez et al., 1997). This is apparent on sensor level but the problem also persists on source level (Schoffelen and Gross, 2009). To address this problem, usually referred to as “artifact of volume conduction,” and to remove or at least attenuate this artifact, several modifications of coupling measures were suggested exploiting the fact that the mixing is essentially instantaneous.

These two questions, what kind of coupling are we interested in and how do we remove artifacts of volume conduction, led to a large variety of coupling measures. Assuming Gaussian distributed data, all non-linear measures must be functions of coherency, and the main content of this paper is the derivation of these functions. This allows to calculate a non-linear coupling measure with a linear model, i.e., we can calculate from empirical data complex coherency and use the respective function as a prediction for the non-linear measure. The difference of the two is then a measure of non-Gaussianity, and it has the potential to detect new phenomena which could otherwise be masked by linear effects.

We here analyze six of these measures plus some variations, namely four non-linear measures of phase-phase coupling, the weighted phase lag index (Vinck et al., 2011), the phase lag index (Stam et al., 2007), the phase locking value, analyzing both the absolute value (Lachaux et al., 1999) and the imaginary part of it (Sadaghiani et al., 2012), and two measures of amplitude-amplitude coupling, one without correction for artifacts of volume conduction (Mehrkanoon et al., 2014) and one with correction for artifacts of volume conduction (Brookes et al., 2012; Hipp et al., 2012).

This paper is organized as follows. We first present background information on linear methods, i.e., coherency and basic functions of it, in section 2.1. In section 2.2 we present the procedure to find or verify mathematical relations numerically. The main part of this paper are sections 2.3 and 2.4 where we present all theoretical findings for phase-phase coupling and amplitude-amplitude coupling, respectively. In section 2.5 we make an excursion to discuss differences between Fourier and Hilbert approach to estimate coupling measures in the frequency domain. We finally present results for simulations, resting state and event related EEG data in sections 3.1, 3.2, and 3.3, respectively. A conclusion is presented in section 4. We tried to keep the main body of the paper as simple as possible, and we therefore moved all mathematical derivations, which are technically quite involved, to an Supplementary Material.

2. Theory

2.1. Background on Linear Coupling Measures

A standard approach to estimate linear relations between two electrophysiological recordings, which can be signals at sensors or estimated sources, as a function of frequency is coherency (Nunez et al., 1997). Typically, data are divided into segments, and for each segment the data are windowed, e.g., using a Hanning window and the Fourier transformations are calculated. Alternative approaches using wavelets or the Hilbert transformation of filtered data are argued to be formally equivalent (Bruns, 2004).

The results are complex numbers z_i(f, k) for the recordings at sensor i at frequency f and segment k. In the following, we will drop the frequency as argument with the implicit understanding that the analysis is done for some given frequency, and we will also omit the segment index k with the implicit understanding that expected values, denoted as < · >, are estimated for empirical data by averaging over k.

For linear and stationary dynamical systems the cross-spectrum contains complete statistical information about the system. It is defined as

$\begin{array}{ll}{S}_{ij}=<z_i{z}_j^{*}>& (1)\end{array}$

where * denotes complex conjugate. Regardless of the details of how data are defined in the frequency domain, for linear and stationary dynamical systems they are always a linear superposition of Gaussian distributed data, are hence themselves Gaussian distributed in the complex domain. Due to stationarity the distribution can only depend on phase differences and not on the phases directly. This distribution is circular Gaussian defined as (Aydore et al., 2013)

$\begin{array}{ll}p(z)=\frac{1}{\text{det}(S){\pi }^2}\text{exp}\left(-z^{†}{S}^{-1}z\right)& (2)\end{array}$

where for a fixed segment index z is the column vector with z_i being its ith element, and † denotes conjugate transpose. The distribution will be used below for all analytic relations between linear and non-linear relations.

The diagonal elements of S are the power values, and the complex coherency C_ij is calculated as

$\begin{array}{ll}{C}_{ij}=\frac{S_{ij}}{{\left(S_{ii}{S}_{jj}\right)}^{1/2}}& (3)\end{array}$

Coherency, like all other measures considered in this paper, can be calculated pairwise. To study relations between different coupling measures it is sufficient to consider only two recordings. For ease of notation, we will therefore omit the sensor index and define coherency c as

$\begin{array}{ll}c=C_{12}=C_{21}^{*}& (4)\end{array}$

Coherency is a complex number. Its absolute value, usually called coherence, is a measure of the strength of the coupling, while its phase is a measure of the average time delay between the peaks of the oscillations. Coherency is a measure of phase-phase coupling, which, however, also depends on amplitude variations because segments of high amplitudes are weighted higher than those with lower amplitude.

The estimation of coupling using EEG and MEG sensor data, and also using respective source estimates, is prone to artifacts of volume conduction. This means that the recordings are mixtures of the true brain activities, and an estimated coupling is likely to be caused by this mixing rather than true coupling between different neuronal sites. To address this problem it was suggested to use the imaginary part of coherency, usually called “imaginary coherence”

$\begin{array}{ll}{c}_I={J}(c)& (5)\end{array}$

where ${J}$(c) denotes imaginary part of c. It can be shown that c_I, also denoted as ImCoh, vanishes for an infinite number of data if all brain sources are independent provided that the quasi-static approximation of the forward model is valid, i.e., the mapping of sources to sensors is instantaneous (Nolte et al., 2004). It should be emphasized that for interacting sources the value of c_I, if non-vanishing, depends on how sources are mapped into sensors.

An important measure for our analysis is lagged coherence (Pascual-Marqui, 2007; Pascual-Marqui et al., 2011), which was proposed both as a signed and unsigned version. We here use the signed version

$\begin{array}{ll}{\tilde{c}}_I=\frac{{J}(c)}{\sqrt{1-{R}{(c)}^2}}& (6)\end{array}$

where ${R}$(c) denotes real part of c.

Lagged coherence, also denoted as LagCoh, has the property that, in contrast to imaginary coherence, its value, apart from sign, does not depend on how the sources are mapped into sources provided that there are only two sources. Also, the possible sign flip is only present in the signed version of lagged coherence and obviously disappears if the absolute value (or square) of ${\tilde{c}}_I$ is taken. We emphasize that this property of being invariant to mixing does not hold for more than two sources. For a proof and simulations of this for two or more sources we refer to Ewald et al. (2012). Therefore, the practical value of this might be limited, because for typical EEG or MEG measurements more than two sources (including all noise sources) are mapped into sensors. In spite of these practical limitations of its interpretation, this quantity seems to play an important and almost universal role for the relation between linear and non-linear coupling measures as we will see below.

2.2. General Remarks on Numerical Evaluation of Coupling Measures

The purpose of simulations is usually to demonstrate the performance of a method under realistic conditions, which we will do below in section 3.1. However, our purpose is to validate mathematical relations with very high accuracy which we would like to illustrate together with theory. To verify analytical results we simulate Gaussian distributed pairs of complex numbers with a random cross-spectral matrix. For each cross-spectral matrix we use 10⁷ pairs of complex numbers. Such a simulation is not meant to represent a realistic measurement but merely to be a tool to find or verify relations with extremely high accuracy. The cross-spectra are constructed as follows. Let

$\begin{array}{ll}x=\left(\begin{array}{c}{x}_1\\ x_2\end{array}\right)=\left(\begin{array}{c}{\eta }_{11}+i{\eta }_{12}\\ {\eta }_{21}+i{\eta }_{22}\end{array}\right)& (7)\end{array}$

where η_nm are independent Gaussian distributed real numbers with zero mean and unit standard deviation and i denotes imaginary unit. Then these numbers are mixed using a random complex mixing matrix A

$\begin{array}{ll}A=A_R+i{A}_I& (8)\end{array}$

where all elements of A_R and A_I are independent Gaussian distributed numbers of zero mean and unit variance. For each mixing matrix A we simulate 10⁷ realizations of observations z as

$\begin{array}{ll}\left(\begin{array}{c}{z}_1\\ z_2\end{array}\right)=z=Ax& (9)\end{array}$

The cross-spectrum of z is then given by

$\begin{array}{ll}S=<z{z}^{†}>=2A{A}^{†}& (10)\end{array}$

Note that in the above formula a factor 2 arises because the unit variances of real and imaginary parts of x add up. Also, S is in general complex because the mixing matrix A is complex which should not be confused with real valued mixing like a mixing artifact occurring in EEG and MEG measurements.

All coupling measures to be analyzed are constructed from expected values of the general form < g(z₁, z₂) >, with the functions g varying across measures, and where those expected values are estimated as averages over all realizations of z.

2.3. Phase-Phase Coupling

2.3.1. wPLI and PLI

Like coherency, the weighted phase lag index (wPLI) is a measure of phase-phase coupling, with averages of phase differences weighted by the amplitudes of individual segments. It is defined here as

$\begin{array}{ll}wPLI=\frac{<{J}(z_1{z}_2^{*})>}{<|{J}(z_1{z}_2^{*})|>}& (11)\end{array}$

Our definition differs slightly from the original definition where in the end the absolute value is taken. Like for lagged coherence, we prefer to keep the sign because the sign might contain relevant information and in general it also simplifies statistics because the absolute value introduces a bias toward positive values. Like lagged coherence, wPLI is invariant to mixing of two sources. This is strictly true only if the absolute value is taken, but the sign may flip using our definition.

We started our tour across a series of coupling measures with wPLI because the relationship to coherency for Gaussian distributed data is known already (Ewald et al., 2012), it has a simple closed form solution, and it illustrates the typical aspects of such relations in general. The relation reads

$\begin{array}{ll}wPLI=\frac{2{\tilde{c}}_I}{1+{\tilde{c}}_I^2}& (12)\end{array}$

where ${\tilde{c}}_I$ is the lagged coherence defined in Equation (6). We consider ${\tilde{c}}_I$ as equivalent to wPLI for Gaussian distributed data in the sense that the latter can be calculated from the former. This is exactly true only for an infinite number of data. For a finite number of data, there are, of course, statistical variations.

We emphasize two important points. First, the equivalence is not trivial, e.g., wPLI is equivalent to ${\tilde{c}}_I$ but not to imaginary coherence c_I. Second, equivalent does not mean identical. The functional relation is simple here, but below we will also see other examples where we observe equivalence clearly from numerical evaluations but the precise functional relation is unclear to us. Numerical results to illustrate these findings are presented in the upper row of Figure 1.

FIGURE 1

Figure 1. Upper row: wPLI as a function of linear measures as indicated for a 500 simulated Gaussian distributed random data sets each consisting of 10⁷ realizations. Lower panels: the same for PLI.

The phase lag index (PLI) is defined as (Stam et al., 2007)

$\begin{array}{ll}PLI=<\text{sign}({J}(z_1{z}_2^{*}))>& (13)\end{array}$

In spite of its name, it is conceptually only loosely related to wPLI. The idea of PLI, using only the sign of the phase differences, is that the resulting coupling should be made less dependent on the actual phase difference. In practice, however, this is hardly the case because for small phase differences a sign flip is more likely than for large phase differences. In contrast to wPLI, PLI is a pure phase measure: any dependence on amplitude was removed in its definition.

As was shown by Pascual-Marqui et al. (2018), PLI, like wPLI, is invariant to mixing the signals within two sensors (or estimated sources) implying that it is invariant of the mixing of sources into sensors provided that there are not more than two sources. This can be used to analytically derive the relation between PLI and linear connectivity measures

$\begin{array}{ll}PLI={\tilde{c}}_I& (14)\end{array}$

The derivation of this is presented in the Supplementary Material. Numerical results to illustrate these findings are presented in the lower row of Figure 1. We emphasize that such an identity is (apparently) true for an infinite number of Gaussian distributed data. For a finite data size, different measures have different statistical properties. The analysis of that is beyond the scope of this paper.

2.3.2. PLV

The phase locking value (PLV) is a classical measure of phase-phase coupling (Lachaux et al., 1999). We here slightly deviate from the original definition, where the absolute value is taken in the end, and define PLV as a complex number:

$\begin{array}{ll}PLV=<\text{exp}(i({\Phi }_1-{\Phi }_2))>=<\frac{z_1{z}_2^{*}}{|z_1||z_2|}>& (15)\end{array}$

with z_k = r_k exp(iΦ_k). Like PLI and in contrast to wPLI it only depends on phase differences and all amplitude variations are ignored.

In Sadaghiani et al. (2012), Bruna et al. (2018), and Palva et al. (2018) it was suggested to use the imaginary part of the complex definition of PLV, referred to as ImPLV in corresponding figures, to construct a measure robust to artifacts of volume conduction. We therefore prefer to keep the complex formulation and present the theory as a whole.

In Aydore et al. (2013) the relation between PLV and coherency was found in terms of hypergeometric functions using Mathematica to solve some of the integrals. In the Supplementary Material we present an explicit derivation which we found numerically to be equivalent to the solution by Aydore et al. (2013). We get the following relation for Gaussian distributed data:

$\begin{array}{ll}PLV=cf(|c|)& (16)\end{array}$

introducing a “scaling function” f which only depends on the absolute value of coherency. We could calculate f analytically only as a series expansion, but not in closed form:

$\begin{array}{ll}f(|c|)=\pi (1-|c{|}^2){\displaystyle \sum _{k=1}^{\infty }}k{\left(\frac{(2k-1)!!}{k!2^k}\right)}^2|c{|}^{2k-2}& (17)\end{array}$

This expansion converges poorly if |c| is close to 1. We refer the reader to the Supplementary Material for an alternative (and less compact) formulation with better convergence properties. There we also give recommendations how to evaluate the function numerically.

We found that f is approximately linear as a function of $\sqrt{1-|c{|}^2}$, and f can be approximated very well by a function $\tilde{f}$ using such a linear function with exact values at the boundaries, namely f (0) = π/4 and f (1) = 1, leading to

$\begin{array}{ll}f(|c|)\approx \tilde{f}(|c|)=1-\left(1-\frac{\pi }{4}\right)\sqrt{1-|c{|}^2}& (18)\end{array}$

Using this approximate function, errors for PLV are smaller than 0.012 which we consider as negligible for practical applications. The functions f and $\tilde{f}$ are shown in Figure 2. For all further analysis we will use the approximate scaling function $\tilde{f}$.

FIGURE 2

Figure 2. Scaling functions f and $\tilde{f}$ as a function of coherence. For the calculation of f(|c|) we used the first 500 terms of the expansion given in Equation (78) of the Appendix.

Numerical results showing the absolute value of PLV as a function of coherence and its scaled version, and also the imaginary part of PLV as a function of imaginary coherence and its scaled version are shown in Figure 3. We observe nearly exact identities for the two scaled versions.

FIGURE 3

Figure 3. Upper row: the absolute value of PLV as a function of coherence (upper left) and as a function of scaled coherence (upper right) for a 500 simulated Gaussian distributed random data sets each consisting of 10⁷ realizations. Lower row: ImPLV as a function of ImCoh (left) and scaled ImCoh (right) for the same random data sets.

2.4. Amplitude-Amplitude Coupling

For amplitude-amplitude coupling relative phase differences are, in its original version, ignored, and the question is whether the amplitudes of two oscillation are functionally related. In principle, the two oscillations can have different frequencies. In such a case we observe cross-frequency coupling, which is always inconsistent with linear dynamical systems. Therefore, we here only consider amplitude-amplitude coupling within a specific frequency.

A functional relation between amplitudes can be measured by a correlation, but details can vary depending on what exactly is correlated. We here consider three different versions: (a) the correlation between powers (the square of the amplitudes), denoted as $corr(|z_1{|}^2,|z_2{|}^2)$ (Mehrkanoon et al., 2014; Soto et al., 2016), (b) the correlation of amplitudes, i.e., corr(|z₁|, |z₂|) (Brookes et al., 2012), and (c) the correlation of the logarithms of the amplitudes, i.e., corr(log|z₁|, log|z₂|) (Hipp et al., 2012). Note, that the latter is equivalent to the correlation of the logarithm of the powers. We here use the term power envelope correlation (PEC) for all these three variants.

We have an analytic solution for Gaussian distributed data only for the first version, which reads

$\begin{array}{ll}corr(|z_1{|}^2,|z_2{|}^2)=|c{|}^2& (19)\end{array}$

where |c| is the coherence of the two signals. A proof is given in the Supplementary Material. A specific consequence is that powers cannot be negatively correlated for linear dynamical systems. Note, that means are subtracted for the calculation of a correlation, and this non-negativity is not a trivial consequence of the positivity of the powers.

Like PLV and coherence, PEC is prone to artifacts of volume conduction. To address this problem it was suggested to replace z₂ by z₂ − αz₁ where the real valued coefficient α is found from fitting αz₁ to z₂ (Brookes et al., 2012, 2014). In the language of the original time series, this means that only that part of the second time series is evaluated which is orthogonal to the first. Without loss of generality, z₁ and z₂ can be normalized such that

$\begin{array}{ll}<|z_1{|}^2>=<|z_2{|}^2>=1& (20)\end{array}$

and then it is straightforward to show that

$\begin{array}{ll}\alpha =c_R& (21)\end{array}$

i.e., α is equal to the real part of coherency.

A different approach was proposed by Hipp et al. (2012) where it was suggested to replace |z₂| by $|{J}(z_2{z}_1^{*})/|z_1||$. The essential difference between these two approaches is that in the first approach α is found globally, i.e., it is the same coefficient for all segments, whereas the second approach is equivalent to fitting a coefficient α separately for each segment. We refer here to the latter approach as a local orthogonalization, and for all of these variants we use the generic term “orthogonalized power envelope correlation” (OPEC).

Similar to PEC, also for OPEC the correlation can refer to power, amplitude or logarithm of the amplitude. We have an analytic solution only for OPEC using power and for global orthogonalization. For normalized signals it reads

$\begin{array}{ll}corr(|z_1{|}^2,|z_2-c_R{z}_1{|}^2)={\tilde{c}}_I^2& (22)\end{array}$

with ${\tilde{c}}_I$ being lagged coherence given in Equation (6). The derivation of this is given in the Supplementary Material.

Numerical results for all nine combinations (three for PEC, three for OPEC with global orthogonalization, and three for OPEC with local orthogonalization) are shown in Figure 4. Remarkably, for all six OPEC variants, the coupling is apparently a unique function of lagged coherence (and then not of imaginary coherence). We also found that correlations of amplitudes are almost identical to those of powers, while the logarithmic transformation has a larger impact.

FIGURE 4

Figure 4. Upper row: three different versions of PEC as a function of the square of coherence for a 500 simulated Gaussian distributed random data sets each consisting of 10⁷ realizations. Middle and lower row: different versions of OPEC as indicated as a function of squared lagged coherence for the same data sets.

2.5. Excursion: Hilbert Transformation vs. Fourier Transformation

In this paper we used the Fourier transformation of windowed time segments to estimate in particular amplitudes as a function of time. A different common practice is to filter the data and then to apply a Hilbert transform to get a complex representation with its absolute value defining the amplitude of the signal. If the original data are Gaussian distributed and stationary, both approaches, being linear and stationary operations, lead to circular Gaussian distributions as discussed, and the theory applies to either approach. On the other hand, explicit results of the two approaches are not identical, and we therefore make a slight excursion to discuss the differences. In Bruns (2004) it is argued that the two approaches are formally equivalent, and any differences are due to typical choices for the windows or filters in either approach. Even though we disagree with this claim as an accurate theorem, we argue that differences between the two approaches are negligible for practical purposes.

The theory for this can be found in Zandvoort and Nolte (2020) and is presented here in more detail. In the Fourier approach the complex signal at frequency f₀ and time t, the time of the center of a corresponding segment, is calculated as

$\begin{array}{ll}y(t,f_0)=(x*h(·,f_0))(t)& (23)\end{array}$

where * denotes convolution and h(t, f₀) is a wavelet. A remark on notation: h(·, f₀) refers to h as the function of time for a specific frequency f₀ to be distinguished from the function value h(t, f₀).

To simplify the notation, we consider an odd number of discrete time points running from −N to N, and the Fourier transforms, to be used below, are defined for discrete frequencies f also running from −N to N. We also assume that the data have zero mean. With such a convention the wavelet is defined as

$\begin{array}{ll}h(t,f_0)=h_0(t)\text{exp}(-i2\pi f_{0}t/(2N+1))& (24)\end{array}$

with h₀(t) being, e.g., the Hanning window. Now, let $\widehat{h}$₀(f) be the Fourier transform of the window, which is substantially different from zero only for small frequencies. Then the Fourier transform of the wavelet reads

$\begin{array}{ll}\widehat{h}(f,f_0)={\widehat{h}}_0(f-f_0)& (25)\end{array}$

with the understanding that we have periodic boundary conditions in the Fourier domain such that $\widehat{h}$₀(f − f₀) is defined for any value of f − f₀. The convolution in Equation (23) is a product in the Fourier domain, and hence

$\begin{array}{ll}y(t,f_0)=\frac{1}{2N+1}{\displaystyle \sum _f}{\widehat{h}}_0(f-f_0)\widehat{x}(f)\text{exp}(i2\pi ft/(2N+1))& (26)\end{array}$

where $\widehat{x}(f)$ is the Fourier transform of x(t). The crucial point now is that, while $\widehat{h}$₀(f) does not vanish for negative frequencies, $\widehat{h}$₀(f − f₀) is negligible for negative frequencies provided that the high frequency f₀ is remote from zero and the Nyquist frequency relative to the width of the wavelet in the frequency domain. In this approximation y(t, f₀) is expressed as an inverse Fourier transform with only positive frequencies as is the case for the Hilbert transform.

In the approach using a Hilbert transform, data are first filtered using a real valued filter. We here focus on conceptual differences between the Fourier and Hilbert approach rather than differences introduced by different filter choices. We hence choose for the Hilbert approach as filter, F(t, f₀), the real part of h(t, f₀)

$\begin{array}{l}F(t,f_0)=\frac{1}{2}{h}_0(t)(\text{exp}(-i2\pi f_{0}t/(2N+1))\\ +\text{exp}(i2\pi f_{0}t/(2N+1)))& (27)\end{array}$

and the filtered data, x_F(t, f₀) = (F(·, f₀) ∗ x(·))(t), can be written as an inverse Fourier transform

$\begin{array}{l}{x}_F(t,f_0)=\frac{1}{2(2N+1)}{\displaystyle \sum _f({\widehat{h}}_0(}f-f_0)\\ +{\widehat{h}}_0(f+f_0))\widehat{x}(f)\text{exp}(i2\pi ft/(2N+1))& (28)\end{array}$

We get the (complex) Hilbert transform of x_F(t, f₀) by summing only over positive frequencies and multiplying the Fourier coefficients with a factor 2:

$\begin{array}{l}{x}_H(t,f_0)=\frac{1}{2N+1}{\displaystyle \sum _{f>0}({\widehat{h}}_0(}f-f_0)\\ +{\widehat{h}}_0(f+f_0))\widehat{x}(f)\text{exp}(i2\pi ft/(2N+1))& (29)\end{array}$

We observe, that the Hilbert transform x_H(t, f₀) is identical to the Fourier formulation y(t, f₀) if the frequency f₀ is remote from zero and the Nyquist frequency relative to the width of the filter, i.e., $\widehat{h}$₀(f − f₀) = 0 for f < 0 and $\widehat{h}$₀(f + f₀) = 0 for f > 0. In practice, this is usually the case to very good approximation.

To illustrate that, we have analyzed empirical resting state EEG data (see below for a description of these data), using both the Hilbert and the Fourier approach to estimate the envelopes. One channel of EEG data containing a strong alpha rhythm at 11 Hz was narrow band filtered at 11 Hz with a FIR filter and its envelope was calculated from the Hilbert transform of the filtered signal. For the Fourier approach we used segments of 1 s duration, windowed the data with a Hanning window, and calculated the Fourier transform for each segment. For this comparison, segments had a maximal overlap, i.e., neighboring segments differed by 1 sample. The Fourier and Hilbert approach can differ in absolute scale depending on the choice of windows and the normalization convention of the Fourier transformation. For simplicity we here have normalized the envelope calculated from the Fourier approach to adjust to the Hilbert envelope. Such global scale factors are irrelevant for all coupling measures studied in this paper. In the left panel of Figure 5 we see for this example that the envelopes for the two approaches are almost identical.

FIGURE 5

Figure 5. (Left) Estimated amplitudes for Hilbert and Fourier approach for an illustrative example. (Middle) Illustration for the same data showing that real and imaginary parts of the Fourier approach are almost perfectly linked by the Hilbert transformation. Non-negligible deviations can be observed in the first and last second of the data, which were taken out from the illustration. (Right) Deviation Δ(f), defined in Equation (30), as a function of the center frequency f of a filter with 1 Hz width applied on white noise. The deviation measures to what accuracy real and imaginary parts of the complex Fourier signals are linked by the Hilbert transformation.

The Hilbert transform links real and imaginary parts of the complex signal. If H denotes the Hilbert transform, then ${J}$(H(${R}$(x))) = ${J}$(x), if x(t) is the complex signal calculated as the Hilbert transform of some real valued signal. In the Fourier approach we didn't calculate the complex signal with a Hilbert transform, and we check whether real and imaginary parts of the complex Fourier representation are linked by a Hilbert transformation. In the middle column we show ${J}$(H(${R}$(y))) vs. ${J}$(y) with y(t) being the complex signal found with the Fourier approach. We observe that those functions are almost identical. Differences can be found at the boundaries of the data set, here the first and last second, which were taken out from the illustration.

As a more general analysis, we calculated the Fourier representation as a function of time for white noise for all frequencies, again using a Hanning window of 1 s duration for all possible temporal positions of the window. To avoid boundary effects, we now strictly imposed periodic boundary conditions such that, e.g., the segment for the Fourier coefficient of the first time point includes data from the end of the whole data set. Let y(t, f) be the Fourier coefficient for frequency f and for the segment with center at time t, then we calculated a deviation Δ(f) as

$\begin{array}{ll}\Delta (f)=\frac{\Vert {J}(y(·,f))-{J}(H({R}(y(·,f)))){\Vert }_F}{\Vert {J}(y(·,f)){\Vert }_F}& (30)\end{array}$

with ‖·‖_F being the Frobenius norm of the corresponding time series. In the right panel of Figure 5 we show Δ(f) as a function of frequency. We generally observe very small deviations apart from moderate deviations for cases where the center frequency of the filter is close to zero or the Nyquist frequency.

3. Results

3.1. Simulations

We simulated both linear and non-linear dynamical systems. Linear dynamics was modeled with an autoregressive model as

$\begin{array}{ll}x(t)={\displaystyle \sum _{\tau =1}^P}A(\tau )x(t-\tau )+\eta (t)& (31)\end{array}$

where for each discrete time point x(t) is a 2 × 1-vector for two channels, the order P was set to 5 and η was independent and white Gaussian distributed noise with zero mean and unit variance. For each data set, all coefficients of the 2 × 2 AR-matrices were set to fixed numbers randomly set from a Gaussian distribution with zero mean and a standard deviation of 0.25. Also, for each set of AR-matrices we checked analytically whether the system was stable and excluded sets with unstable dynamics. The physical units can be arbitrarily set, and we defined the sampling rate to be 100 Hz. Also, to concentrate on systematic effects, we simulated rather long data sets of 30 min. For a total of 100 data sets we analyzed coupling always at 10 Hz. Here and in the following, data were divided onto segments of 1 s duration with 50% overlap. The data in each segment were Hanning-windowed, and Fourier coefficients at 10 Hz were used as input for all coupling estimates.

As a non-linear dynamical system we used the Kuramoto-model chosen here with stochastic input, random coupling parameters and random delays. Specifically, the dynamics was defined for discrete time t for two channels i as

$\begin{array}{l}\theta (t+1,i)=\theta (t,i)+dt(\omega +\sigma \eta (t,i)+{\displaystyle \sum _{j=1}^{2}a}(i,j)\text{ sin}(\theta (t,i)\\ -\theta (t-\tau (i,j),j)))& (32)\end{array}$

with the settings dt = 1/100 corresponding to a sampling rate of 100 Hz, ω = 20π resulting in 10 Hz oscillations, and σ = 5. The delays τ(i, j) were set at (time-independent) random integer numbers corresponding to delays up to 100 ms. The coupling parameters were set to a(i, i) = 0, and a(i, j) for i ≠ j was set to a random number from a Gaussian distribution with zero mean and standard deviation 0.5. Results are not crucial with respect to details of the simulations, but our choice resulted into a roughly even coverage of all possible outcomes for the linear coupling measures. From the phases the time series were finally calculated as

$\begin{array}{ll}x(t,i)=\text{sin}(\theta (t,i))& (33)\end{array}$

for time t and channel i. We note that even though we did not introduce an explicit time dependent amplitude, the amplitude estimate varied in time because the phases were not linear functions of time leading to variations of the instantaneous frequency and then to variations of the amplitude estimate at a specific frequency. As for the linear dynamics we simulated 100 data sets with 30 min of data each.

Results for both the linear and non-linear dynamics are shown in Figure 6 where we display the non-linear coupling value for six different coupling measures as a function of the corresponding value calculated from the linear coupling values using the derived equations for Gaussian distributed data. For the linear dynamical systems we observe nearly perfect correspondence between the two. For the non-linear systems and for phase-phase coupling we observe only minor but clearly systematic deviations. For the non-linear systems, amplitudes are almost not coupled at all and a prediction of that coupling assuming Gaussian distributed data fails completely.

FIGURE 6

Figure 6. Illustration of non-linear coupling measures as a function of their corresponding model values calculated from complex coherency and assuming Gaussian distributed data for simulated data. Blue dots refer to linear dynamical systems, and red dots to the Kuramoto models each for a different specification of parameters.

3.2. Resting State EEG

We analyzed cleaned resting state EEG data measured with eyes closed for 10 subjects publicly available at http://clopinet.com/causality/data/nolte/. The data consist of around 10 min recordings in 19 channels with mathematically linked ears reference. The data are used here such that our results can be reproduced. Our complete code for the analysis is available upon request. The data are a subset of data for 88 subjects, which are described in more detail in Nolte et al. (2008). Only this subset is publicly available.

First of all, for these data sets we analyzed how well non-linear coupling matrices can be explained by the respective linear models. Let D(f, k) be a connectivity matrix for all pairs of sensors calculated with a specific measure for frequency f and subject k, and let D_M(f, k) be the corresponding model connectivity matrix calculated from the coherency matrix. In Figures 7, 8 we show two illustrative examples, calculated from the first subject at frequencies 10 and 15 Hz, respectively. For each pair of sensors we show the actual non-linear measure and the result of the corresponding model. Resting state with eyes closed is known to show a strong 10 Hz rhythm (the alpha rhythm) consisting of interacting sources with substantial phase delays. Such interactions are also observable at other frequencies, but typically to a much lesser extent.

FIGURE 7

Figure 7. Illustration of non-linear coupling measures as a function of their corresponding model values calculated from complex coherency and assuming Gaussian distributed data. Each dot refers to one pair of sensors for one subject with coupling values calculated at 10 Hz. For PEC and OPEC, the non-linear measures can also come out negative, but that's very rare, values are typically only slightly below zero, and such a case did not occur in this example.

FIGURE 8

Figure 8. Same as Figure 7 but now for coupling measures calculated at 15 Hz. We note again that PEC and OPEC can be negative, but such a case did not occur in this example.

We consider the following observations as most relevant. The wPLI can be explained by the model almost perfectly, regardless of whether there are substantial interactions (with delay) or not, amplitude-amplitude coupling (i.e., PEC) is systematically larger than the model prediction, and amplitude-amplitude coupling with attenuation of artifacts of volume conduction (i.e., OPEC) cannot be explained by the model at all. This could be a problem when interpreting OPEC as a coupling measure robust to artifacts of volume conduction as will be discussed in the conclusion in section 4.

In the following we present a systematic analysis including results for all subjects and frequencies. For K subjects we define an average model error for each frequency as

$\begin{array}{ll}{E}_M(f)=\frac{1}{K}{\displaystyle \sum _{k=1}^K}\frac{\Vert D(f,k)-D_M(f,k){\Vert }_F}{\Vert D(f,k){\Vert }_F}& (34)\end{array}$

where ‖·‖_F denotes Frobenius norm.

Non-vanishing model errors can have two causes: (a) the data are non-Gaussian distributed, and (b) they are caused by statistical fluctuations. To assess the magnitude of the statistical fluctuation we also calculated a statistical error. For this we replaced the model connectivity by a connectivity matrix calculated with the non-linear measure with resampled data using the bootstrap procedure where we constructed new data of the same size by randomly picking segments of the original data with repetition. The bootstrap resampling generates data from approximately the same distribution as the original data without destroying any coupling and then serves to estimate statistical variations in a generic way. For the k.th subject and frequency f we calculated N = 20 such connectivity matrices denoted as D_S(f, k, n) for n = 1..N, and a statistical error was calculated as

$\begin{array}{ll}{E}_S(f)=\frac{1}{KN}{\displaystyle \sum _{k,n}}\frac{\Vert D(f,k)-D_S(f,k,n){\Vert }_F}{\Vert D(f,k){\Vert }_F}& (35)\end{array}$

If this statistical error is larger than or of similar size as the model error, this model error should not be interpreted as an indicator of non-Gaussianity of the data. Results for six different non-linear coupling measures are shown in Figure 9. In addition, we also calculated the statistical errors for the linear measures. We observe that generally for phase-phase coupling all methods with correction for artifacts of conduction are similar: results are statistically unstable for frequencies outside the alpha band. For the frequencies outside the alpha band the models are typically poor which is probably not surprising as coupling measures robust to artifacts of volume conduction are relatively weak and can hardly be reproduced with different methods. An exception is wPLI which can always be explained very well with the linear model, indicating that wPLI depends very little on non-linear properties of the data. The absolute value of PLV can typically be explained very well with the linear model, but to lesser extent in the alpha band. This is the only phenomenon for phase-phase coupling, where we can clearly observe deviations of the model and actual coupling larger than the statistical error.

FIGURE 9

Figure 9. Upper and middle row: relative model errors (full lines) averaged over 10 subjects for six different non-linear coupling measures each of them modeled by the corresponding function of coherency assuming Gaussian distributed data. Estimated statistical errors are shown as dashed lines. Lower row: statistical errors for three different linear measures. We emphasize that scales for the y-axis vary across subplots.

Systematic deviations of the model predictions larger than statistical errors can be observed for PEC and OPEC for all frequencies. While for PEC the model (which is the square of coherence) explains around 80% of the observation (model error is around 0.2), the model error for OPEC is nearly 100% which is not surprising when inspecting the examples shown in Figures 7, 8.

3.3. Reaction Task With and Without Conflict

We now present results for an event related paradigm. The purpose here is merely an illustration to show that we can get non-trivial results also in such a case. Details of the experiment can be found in Li et al. (2015), and we here only give a short description. In this paradigm, the stimulus was an upward or downward arrow presented at one of four possible locations of the screen: top left, top right, bottom left, and bottom right. Participants were asked to respond to the direction of the arrow as soon as possible by pressing the “F” key or the “J” key on a keyboard, while ignoring the location of the arrow, which was either congruent or incongruent with the direction of the arrow. The mapping of arrow direction and response key was counterbalanced between participants. EEG was measured in 62 channels (plus 2 mastoids, which were not included in the connectivity analysis) and referenced to the mathematically linked mastoids. In total, we analyzed 33 subjects with an average of 228 trials per condition. For each trial, 1,200 ms of data from 200 ms before the stimulus until 1 s after the stimulus were analyzed further. ERPs were subtracted from the raw data such that a connectivity analysis corresponds to the analysis of fluctuations around the ERPs. Each trial was divided into segments of 200 ms duration with an overlap of 180 ms such that we could calculate connectivity for 51 different time points. Of course, such short segments of 200 ms result in a poor frequency resolution of Δf ≈ 5 Hz, but with a high frequency resolution time dependence cannot be analyzed anymore.

For these data, we only analyzed PEC and calculated the correlations of the squares of the amplitudes. We recall that the model assuming Gaussian distributed data predicts that PEC is the square of coherence. Similar to the analysis of the resting state data we calculated a model error with Equation (34), which now also depends on the time of the segments relative to the stimulus. We observed that results are very similar for conflict and non-conflict trials (not shown), and we therefore present only results where we combined the conditions. The model error is shown in the upper left panel of Figure 10. The model error is relatively large (around 0.2) for the alpha range before the stimulus and at the end of trials. It drops substantially in the center of the trials.

FIGURE 10

Figure 10. Upper row: model errors as a function of time and frequency for PEC using correlation of powers. The stimulus was at t = 0ms. In the left panel we show original model errors, and in the right panel we show model errors with results at baseline (the first 200 ms) subtracted, and with non-significant values set to zero. Lower panels: the same for powers, but here we always normalize power values by the power values at the baseline.

The question here is whether changes in time can be detected significantly. Therefore, we calculated the difference of the model error to the baseline, which we set to be the results at the first time point, i.e., the segments before the stimulus. Significance was tested using a paired permutation test: for each subject and for each time point, results (at the same frequency) were randomly switched between the baseline and the actual time point. We constructed 10,000 such surrogates and the p-value was calculated as the fraction of cases for which the surrogates showed a larger difference of baseline and actual time point than the original data. We analyzed nine frequencies from 5 to 45 Hz. Using the Bonferroni correction for 51 time points and these nine frequencies we considered results as significant if the p-value was lower than 0.05/(51 ∗ 9). In the upper right panel of Figure 10 we show the model error with baseline subtracted and after setting non-significant differences to zero. We observe significant time variation of the model error mainly in the alpha range.

In the lower panels of Figure 10 we show analogous results for power, here always showing the power ratio of the power for each time point and the corresponding power at baseline. While we observe a drop in alpha power (lower left panel), this drop is not significant.

4. Conclusion

In this paper we presented mathematical relations between linear and non-linear measures of brain coupling assuming Gaussian distributed data. All relations were verified in simulations. Let us recall the main theoretical results. We considered four different non-linear measures of phase-phase coupling: wPLI, PLI, and absolute value and imaginary part of PLV. The functional relations could be proven for all of these measures. For amplitude-amplitude coupling we could solve the problem analytically in closed form only if powers (rather than amplitudes or the logarithm of powers) are correlated and if artifacts of volume conduction are corrected for globally (i.e., time independent) or not at all. All other variants could only be analyzed numerically. To our own surprise we found that all considered versions of amplitude-amplitude coupling with correction for mixing artifacts turned out to be functions of lagged coherence. Except for the one case we could solve analytically, we must leave the mathematical proof of this as an open question. In total, with the exception of the imaginary part of PLV, all considered measures with correction for mixing artifacts turned out to be functions of lagged coherence for Gaussian distributed data, and hence lagged coherence plays a quite universal role.

In general, calculating linear measures from the cross-spectrum has a couple of computational advantages. First, if the data are Gaussian distributed, then the cross-spectrum is a maximum likelihood estimator of the parameters of the distribution and then the estimation has the smallest statistical error. Second, when using a linear inverse method (or a quasi-linear method like a beamformer) the cross-spectrum in source space can easily be calculated from the cross-spectrum in sensor space by multiplying the latter with the spatial filter from the left and right. In contrast, when using a non-linear method, the entire raw data need to be mapped into source space. Third, in case of amplitude-amplitude coupling phase information contained in complex coherency gets lost. The question then is whether for empirical data the calculation of non-linear measures in particular in source space is worth the effort and adds any useful information. Empirical electrophysiological data are certainly not Gaussian distributed, even if the deviation of Gaussianity is usually weak. The mathematical relations presented here may serve as a tool to isolate non-linear effects by subtracting from the non-linear measure the corresponding model results assuming a Gaussian distribution. Such a difference cannot be explained by any linear dynamical model. If, e.g., a substantial non-vanishing difference is observed, this could indicate that a linear autoregressive model is inadequate to calculate Granger Causality.

For empirical data we analyzed how well observed non-linear coupling measures could be explained by calculating coherency from these data and then predicting the non-linear coupling measure assuming Gaussian distributed data and using the theoretical relations. For event related data, we only illustrated the procedure and could show that deviations of the model for amplitude-amplitude coupling are in general time dependent. A more complete analysis was given for resting state EEG data. For phase-phase coupling we found that deviations are minor, i.e., essentially within statistical errors, with the exception of the absolute value of PLV at 10 Hz. Deviations were quite substantial for amplitude-amplitude coupling, in particular when a correction for artifacts of volume conduction is included. That measure is known to be robust to artifacts of volume conduction in general only for Gaussian distributed data, as made clear in much detail by the authors themselves (Brookes et al., 2014). To give an independent explanation of this we shortly discuss conceptual relations between OPEC and independent component analysis (ICA). The ICA model assumes that any coupling between sensors is a mixing artifact and ICA attempts to demix to independent components. ICA gives a unique decomposition only when sources are non-Gaussian distributed. OPEC fails for non-Gaussian distributed data and ICA works only for non-Gaussian distributed for the same reason: for non-Gaussian distributed data, removing linear dependencies is not sufficient to get rid of mixing artifacts from coupling measures. This is the problem for OPEC but it is exploited in ICA to get additional information on coupling to finally arrive at a unique decomposition.

These conceptual issues raises the question how to justify the usage of OPEC. If we need to assume that data are Gaussian distributed to reasonable approximation to justify a specific conclusion we might as well calculate the linear coupling measure in the first place. In any event, whether the large deviations from the Gaussian model observed for OPEC are artifacts of volume conduction or correspond to genuine non-linear brain interactions remains an open question which needs to be addressed in the future.

This paper leaves many other questions open. Apart from the lack of mathematical proofs for some cases, the analysis of empirical data was rather coarse. We estimated model errors as averages over all sensor pairs. Interesting effects, e.g., delayed brain interactions at other frequencies than 10 Hz, certainly exist but can easily be masked by such averages if those effects are relatively weak and/or occur only in a few sensor pairs. Also, our analysis was done completely in sensor space, and the question remains open where in the brain we observe large or small deviations from the linear model, and, most importantly, whether those differences can be explained by a reasonable model of brain dynamics. Finally, the question is whether we can observe differences of these deviations for different experimental conditions or brain pathologies. All these questions are beyond the scope of this paper and need to be addressed in the future.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found at: http://clopinet.com/causality/data/nolte/.

Ethics Statement

The studies involving human participants were reviewed and approved by Institutional Review Board of the Institute of Psychology, Chinese Academy of Sciences. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

GN developed the theory. GN and EG-L conducted the simulations. XL and ZL provided the event-related EEG data. All authors contributed to the writing.

Funding

This research was partially funded by the BMBF (161A130), the German Research Foundation (DFG, SFB936/A2/A3/Z3, TRR169/B1/B4, and SPP2041/EN533/15-1), and from the Landesforschungsförderung Hamburg (CROSS, FV25).

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.

Supplementary Material

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

References

Aydore, S., Pantazis, D., and Leahy, R. (2013). A note on the phase locking value and its propertie. Neuroimage 74, 231–244. doi: 10.1016/j.neuroimage.2013.02.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Baccala, L., and Sameshima, K. (2001). Partial directed coherence: a new concept in neural structure determination. Biol. Cybern. 84, 463–474. doi: 10.1007/PL00007990

PubMed Abstract | CrossRef Full Text | Google Scholar

Bressler, S. L., and Seth, A. K. (2011). Wiener-Granger causality: a well established methodology. Neuroimage 58, 323–329. doi: 10.1016/j.neuroimage.2010.02.059

PubMed Abstract | CrossRef Full Text | Google Scholar

Brookes, M. J., O'Neill, G. C., Hall, E. L., Woolrich, M. W., Baker, A., Corner, S. P., et al. (2014). Measuring temporal, spectral and spatial changes in electrophysiological brain network connectivity. Neuroimage 91, 282–299. doi: 10.1016/j.neuroimage.2013.12.066

PubMed Abstract | CrossRef Full Text | Google Scholar

Brookes, M. J., Woolrich, M. W., and Barnes, G. R. (2012). Measuring functional connectivity in MEG: a multivariate approach insensitive to linear source leakage. Neuroimage 63, 910–920. doi: 10.1016/j.neuroimage.2012.03.048

PubMed Abstract | CrossRef Full Text | Google Scholar

Bruna, R., Maestu, F., and Pereda, E. (2018). Phase locking value revisited: teaching new tricks to an old dog. J. Neural Eng. 15:056011. doi: 10.1088/1741-2552/aacfe4

PubMed Abstract | CrossRef Full Text | Google Scholar

Bruns, A. (2004). Fourier-, hilbert- and wavelet-based signal analysis: are they really different approaches? J. Neurosci. Methods 137, 321–332. doi: 10.1016/j.jneumeth.2004.03.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Engel, A. K., Fries, P., and Singer, W. (2001). Dynamic predictions: oscillations and synchrony in top-down processing. Nat. Rev. Neurosci. 2, 704–716. doi: 10.1038/35094565

PubMed Abstract | CrossRef Full Text | Google Scholar

Engel, A. K., Gerloff, C., Hilgetag, C. C., and Nolte, G. (2013). Intrinsic coupling modes: multiscale interactions in ongoing brain activity. Neuron 80, 867–886. doi: 10.1016/j.neuron.2013.09.038

PubMed Abstract | CrossRef Full Text | Google Scholar

Ewald, A., Marzetti, L., Zappasodi, F., Meinecke, F. C., and Nolte, G. (2012). Estimating true brain connectivity from EEG/MEG data invariant to linear and static transformations in sensor space. Neuroimage 60, 476–488. doi: 10.1016/j.neuroimage.2011.11.084

PubMed Abstract | CrossRef Full Text | Google Scholar

Fries, P. (2005). A mechanism for cognitive dynamics: neuronal communication through neuronal coherence. Trends Cogn. Sci. 9, 474–480. doi: 10.1016/j.tics.2005.08.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Fries, P. (2015). Rhythms for cognition: communication through coherence. Neuron 88, 220–235. doi: 10.1016/j.neuron.2015.09.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Hipp, J. F., Hawellek, D. J., Corbetta, M., Siegel, M., and Engel, A. K. (2012). Large-scale cortical correlation structure of spontaneous oscillatory activity. Nat. Neurosci. 15, 884–890. doi: 10.1038/nn.3101

PubMed Abstract | CrossRef Full Text | Google Scholar

Kaminski, M. J., and Blinowska, K. (1991). A new method of the description of the information-flow in the brain structures. Biol. Cybern. 65, 203–210. doi: 10.1007/BF00198091

PubMed Abstract | CrossRef Full Text | Google Scholar

Lachaux, J., Rodriguez, E., Martinerie, J., and Varela, F. (1999). Measuring phase synchrony in brain signals. Hum. Brain Mapp. 8, 194–208. doi: 10.1002/(SICI)1097-0193(1999)8:4<194::AID-HBM4>3.0.CO;2-C

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Q., Wang, K., Nan, W., Zheng, Y., Wu, H., Wang, H., et al. (2015). Electrophysiological dynamics reveal distinct processing of stimulus-stimulus and stimulus-response conflicts. Psychophysiology 52, 562–571. doi: 10.1111/psyp.12382

PubMed Abstract | CrossRef Full Text | Google Scholar

Mehrkanoon, S., Breakspear, M., Britz, J., and Boonstra, T. (2014). Intrinsic coupling modes in source-reconstructed electroencephalography. Brain Connect. 10, 812–825. doi: 10.1089/brain.2014.0280

PubMed Abstract | 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, 2292–2307. doi: 10.1016/j.clinph.2004.04.029

PubMed Abstract | CrossRef Full Text | Google Scholar

Nolte, G., Ziehe, A., Nikulin, V., Schlögl, A., Krämer, N., Brismar, T., et al. (2008). Robustly estimating the flow direction of information in complex physical systems. Phys. Rev. Lett. 100:234101. doi: 10.1103/PhysRevLett.100.234101

PubMed Abstract | CrossRef Full Text | Google Scholar

Nunez, P., Srinivasan, R., Westdorf, A., Wijesinghe, R., Tucker, D., Silberstein, R., et al. (1997). EEG coherency. I. Statistics, reference electrode, volume conduction, Laplacians, cortical imaging, and interpretation at multiple scales. Electroencephalogr. Clin. Neurophysiol. 103, 499–515. doi: 10.1016/S0013-4694(97)00066-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Palva, J. M., Wang, S. H., Palva, S., Zhigalov, A., Monto, S., Brookes, M. J., et al. (2018). Ghost interactions in MEG/EEG source space: a note of caution on inter-areal coupling measures. Neuroimage 173, 632–643. doi: 10.1016/j.neuroimage.2018.02.032

PubMed Abstract | CrossRef Full Text | Google Scholar

Pascual-Marqui, R. D. (2007). Instantaneous and lagged measurements of linear and nonlinear dependence between groups of multivariate time series: frequency decomposition. arXiv 0711.1455.

Google Scholar

Pascual-Marqui, R. D., Faber, P., Kinoshita, T., Kochi, K., Milz, P., Nishida, K., et al. (2018). A comparison of bivariate frequency domain measures of electrophysiological connectivity. bioRxiv 459503. doi: 10.1101/459503

CrossRef Full Text | Google Scholar

Pascual-Marqui, R. D., Lehmann, D., Koukkou, M., Kochi, K., Anderer, P., Saletu, B., et al. (2011). Assessing interactions in the brain with exact low-resolution electromagnetic tomography. Phil. Trans. R. Soc. A 369, 3768–3784. doi: 10.1098/rsta.2011.0081

PubMed Abstract | CrossRef Full Text | Google Scholar

Sadaghiani, S., Scheeringa, R., Lehongre, K., Morillon, B., Giraud, A.-L., D'Esposito, M., et al. (2012). Alpha-band phase synchrony is related to activity in the fronto-parietal adaptive control network. J. Neurosci. 32, 14305–14310. doi: 10.1523/JNEUROSCI.1358-12.2012

PubMed Abstract | CrossRef Full Text | Google Scholar

Schoffelen, J., and Gross, J. (2009). Source connectivity analysis with MEG and EEG. Hum. Brain Mapp. 30, 1857–1865. doi: 10.1002/hbm.20745

PubMed Abstract | CrossRef Full Text | Google Scholar

Soto, J. L. P., Lachaux, J. P., Baillet, S., and Jerbi, K. (2016). A multivariate method for estimating cross-frequency neuronal interactions and correcting linear mixing in MEG data, using canonical correlations. J. Neurosci. Methods 271, 169–181. doi: 10.1016/j.jneumeth.2016.07.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Stam, C., Nolte, G., and Daffertshofer, A. (2007). Phase lag index: assessment of functional connectivity from multi channel EEG and MEG with diminished bias from common sources. Hum. Brain Mapp. 28, 1178–1193. doi: 10.1002/hbm.20346

PubMed Abstract | CrossRef Full Text | Google Scholar

Vinck, M., Oostenveld, R., van Wingerden, M., Battaglia, F., and Pennartz, C. M. A. (2011). An improved index of phase-synchronization for electrophysiological data in the presence of volume-conduction, noise and sample-size bias. Neuroimage 55, 1548–1565. doi: 10.1016/j.neuroimage.2011.01.055

PubMed Abstract | CrossRef Full Text | Google Scholar

Zandvoort, C., and Nolte, G. (2020). Understanding phase-amplitude coupling from bispectral analysis. bioRxiv. doi: 10.1101/2020.03.20.000745

CrossRef Full Text | Google Scholar