Brain connectivity estimators

Classical estimators of connectivity are correlation and coherence. The above measures provide information on the directionality of interactions in terms of delay (correlation) or coherence (phase), however the information does not imply causal interaction. Moreover, it may be ambiguous, since phase is determined modulo 2π. It is also not possible to identify by means of correlation or coherence.

The most frequently used nonlinear estimators of connectivity are mutual information, transfer entropy, generalised synchronisation,[7] the continuity measure,[8] synchronization likelihood,[9] and phase synchronization.[7] Mutual information and transfer entropy rely on the construction of histograms for probability estimates. The continuity measure, generalized synchronisations, and synchronisation likelihood are very similar methods based on phase space reconstruction. Among these measures, only transfer entropy allows for the determination of directionality. Nonlinear measures require long stationary segments of signals, are prone to systematic errors, and above all are very sensitive to noise.[7][8][10] The comparison of nonlinear methods with linear correlation in the presence of noise reveals the poorer performance of non-linear estimators.[8] In [7] the authors conclude that there must be good reason to think that there is non-linearity in the data to apply non-linear methods. In fact it was demonstrated by means of surrogate data test,[11][12] and time series forecasting [13] that nonlinearity in EEG and LFP is the exception rather than the norm. On the other hand, linear methods perform quite well for non-linear signals.[14] Finally, non-linear methods are bivariate (calculated pair-wise), which has serious implication on their performance.

The AR model assumes that X(t)—a sample of data at a time t—can be expressed as a sum of p previous values of the samples from the set of k-signals weighted by model coefficients A plus a random value E(t):

${\displaystyle X(t)=\sum _{j=1}^{p}A(j)X(t-j)+E(t)}$

(1)

The p is called the model order. For a k-channel process X(t) and E(t) are vectors of size k and the coefficients A are k×k-sized matrices. The model order may be determined by means of criteria developed in the framework of information theory and the coefficients of the model are found by means of the minimalization of the residual noise. In the procedure correlation matrix between signals is calculated. By the transformation to the frequency domain we get:

${\displaystyle {\begin{array}{l}E(f)=A(f)X(f)\\X(f)=A^{-1}(f)E(f)=H(f)E(f)\end{array}}}$

(2)

H(f) is a transfer matrix of the system, it contains information about the relationships between signals and their spectral characteristics. H(f) is non-symmetric, so it allows for finding causal dependencies. Model order may be found by means of criteria developed in the framework of information theory,[19] e.g. AIC criterion.

Granger causality index showing the driving of channel x by channel y is defined as the logarithm of the ratio of residual variance for one channel to the residual variance of the two-channel model:[21] GCI_y→x = ln (e/e₁) This definition can be extended to the multichannel system by considering how the inclusion of the given channel changes the residual variance ratios. To quantify directed influence from a channel x_j to x_i for n channel autoregressive process in time domain we consider n and n−1 dimensional MVAR models. First, the model is fitted to whole n-channel system, leading to the residual variance V_i,n(t) = var(E_i,n(t)) for signal x_i. Next, a n−1 dimensional MVAR model is fitted for n−1 channels, excluding channel j, which leads to the residual variance V_i,n−1(t) = var (E_i,n−1(t)). Then Granger causality is defined as:

${\displaystyle \mathrm {GCI} _{j\rightarrow i}(t)=\ln \left({\frac {V_{i,n}(t)}{V_{i,n-1}(t)}}\right)}$

GCI is smaller or equal 1, since the variance of n-dimensional system is lower than the residual variance of a smaller, n−1 dimensional system. GCI(t) estimates causality relations in time domain. For brain signals the spectral characteristics of the signals is of interest, because for a given task the increase of propagation in certain frequency band may be accompanied by the decrease in another frequency band.[22] DTF or PDC are the estimators defined in the frequency domain.

Directed Transfer Function (DTF) was introduced by Kaminski and Blinowska [23] in the form:

${\displaystyle \mathrm {DTF} _{j\rightarrow i}^{2}(f)={\frac {\left|H_{ij}(f)\right|^{2}}{\sum _{m=1}^{k}\left|H_{im}(f)\right|^{2}}}}$

(3)

Where H_ij(f) is an element of a transfer matrix of MVAR model. DTF describes causal influence of channel j on channel i at frequency f. The above equation (3) defines a normalized version of DTF, which takes values from 0 to 1 producing a ratio between the inflow from channel j to channel i to all the inflows to channel i. The non-normalized DTF which is directly related to the coupling strength [24] is defined as:

${\displaystyle \mathrm {NDTF} _{j\rightarrow i}^{2}(f)=\left|H_{ij}(f)\right|^{2}}$

(4)

DTF shows not only direct, but also cascade flows, namely in case of propagation 1→2→3 it shows also propagation 1→3. In order to distinguish direct from indirect flows direct Directed Transfer Function (dDTF) was introduced.[25] The dDTF is defined as a multiplication of a modified DTF by partial coherence. The modification of DTF concerned normalization of the function in such a way as to make the denominator independent of frequency. The dDTF_j→i showing direct propagation from channel j to i is defined as:

${\displaystyle {\begin{array}{l}\mathrm {dDTF} _{j\rightarrow i}^{2}(f)=F_{ij}^{2}(f)C_{ij}^{2}(f)\\F_{ij}^{2}(f)=\displaystyle {\frac {\left|H_{ij}(f)\right|^{2}}{\sum _{f}\sum _{m=1}^{k}\left|H_{im}(f)\right|^{2}}}\end{array}}}$

(5)

Where C_ij(f) is partial coherence. The dDTF_j→i has a nonzero value when both functions F_ij(f) and C_ij(f) are non-zero, in that case there exists a direct causal relation between channels j→i. Distinguishing direct from indirect transmission is essential in case of signals from implanted electrodes, for EEG signals recorded by scalp electrodes it is not really important.[15]

DTF may be used for estimation of propagation in case of point processes e.g. spike trains or for the estimation of causal relations between spike trains and Local Field Potentials.[26]

The partial directed coherence (PDC) was defined by Baccala and Sameshima [27] in the following form:

${\displaystyle P_{ij}(f)={\frac {A_{ij}(f)}{\sqrt {\mathbf {a} _{j}^{*}(f)\mathbf {a} _{j}(f)}}}}$

(6)

In the above equation A_ij(f) is an element of A(f)—a Fourier transform of MVAR model coefficients A(t), where a_j(f) is j-th column of A(f) and the asterisk denotes the transpose and complex conjugate operation. Although it is a function operating in the frequency domain, the dependence of A(f) on the frequency has not a direct correspondence to the power spectrum. From normalization condition it follows that PDC takes values from the interval [0,1]. PDC shows only direct flows between channels. Unlike DTF, PDC is normalized to show a ratio between the outflow from channel j to channel i to all the outflows from the source channel j, so it emphasizes rather the sinks, not the sources. The normalization of PDC affects the detected intensities of flow as was pointed out in.[28] Namely, adding further variables that are influenced by a source variable decreases PDC, although the relationship between source and target processes remains unchanged. In other words: the flow emitted in one direction will be enhanced in comparison to the flows of the same intensity emitted from a given source in several directions.