Evaluation of Phase Locking and Cross Correlation Methods for Estimating the Time Lag between Brain Sites: A Simulation Approach.

INTRODUCTION: Direction and latency of electrical connectivity between different sites of brain explains brain neural functionality. We compared efficiency of cross correlation and phase locking methods in time lag estimation which are based on local field potential (LFP) and LFP-spike signals, respectively.
METHODS: Signals recorded from MT area of a macaque's brain was used in a simulation approach. The first signal was real brain activity and the second was identical to the first one, but with two kinds of delayed and not delayed forms. Time lag between two signals was estimated by cross correlation and phase locking methods.
RESULTS: Both methods estimated the time lags with no errors. Phase locking was not as time efficient as correlation. In addition, phase locking suffered from temporal self bias. DISCUSSION: Correlation was a more efficient method. Phase locking was not considered as a proper method to estimate the time lags between brain sites due to time inefficiency and self bias, the problems which are reported for the first time about this method.

1. Introduction

The local field potentials (LFP) can be recorded by multi-electrodes from various sites of neural tissue of the brain. To study the brain functionality, it is important to discover the direction and velocity of activation of brain sites during a sensory/motor task. Various methods have been proposed in the literature to determine the direction of brain activity, i.e. whether the activity of a site takes place before or after the other, as well as to estimate the time latency after which the activation occurs (indicating propagation velocity). Some methods are limited to direction of recognition rather than estimation of time lag between regions, while others estimate the time lag (by means of which the direction is also extractable). Studies such as (Alonso, and Martinez, 1998; Holdefer, Miller, Chen, & Houk, 2000; Lindsey, Hernandez, Morris, Shannon, 1992; Snider, Kabara, Roig, & Bonds, 1998) estimated the time lag by estimating cross correlation of spike trains. Also some other methods such as covariance were applied based on spike activity (Paz, Bauer, & Pare, 2009; Siapas, Lubenov, & Wilson, 2005). Since it is not possible to record spikes in some cases, such methods are not always useful. Granger causality (Cadotte, Demarse, Mareci, Parekh, Talathi, & Hwang et al., 2010; Gregoriou, Gotts, Zhou, & Desimone, 2009; Popa, Duvarci, Popescu, Lena, & Pare, 2010) and partial directed coherence (PDC) (Astolfi, Cincotti, Mattia, Marciani, Baccala, & de Vico Fallani et al., 2006; Baccala, & Sameshima, 2001; Taxidis, Coomber, Mason, & Owen, 2010; Winterhalder, Schelter, Hesse, Schwab, Leistritz, L., & Klan et al., 2005) can discover the direction by the LFP's. They are rather complicated (Gourevitch, Bouquin-Jeannes, & Faucon, 2006) and sensitive to noise (Taxidis, Coomber, Mason, & Owen, 2010; Winterhalder, Schelter, Hesse, Schwab, Leistritz, L., & Klan et al., 2005). Cross correlation of LFP's (Adhikari, Sigurdsson, Topiwala, & Gordon, 2010) is a simple method to estimate the time lag which proved to have lower sensitivity to noise in comparison with PDC in direction of recognition. Time-delayed mutual information of activation phase (Wilmer, Lussanet, & Lappe, 2012) led to an information-theoretic estimation of lag between brain sites. The main feature of the method is assuming nonlinear relation between the brain sites under the study; nevertheless, it is not as simple as cross correlation method.

As indicated in a study (Onslow, Bogacz, & Jones, 2011), varying data length and noise levels may lead to false positives measuring of phase-amplitude coupling (PAC) by envelope-to-signal correlation (ESC), the modulation index (MI) and cross-frequency coherence (CFC). (Cohen, Bour, Mantione, Figee, Vink, Tijssen, Rootselaar, Munckhof, Schuurman, & Denys, 2011) applied spectral Granger causality analyses to discover the direction of functional connectivity between medial frontal cortex (MFC) and nucleus accumbens and demonstrated the role of the MFC in modulation of nucleus accumbens processing. A time-varying instantaneous frequency detection approach (Zaen, Murray, Meuli, & Vesin, 2013) was proposed to estimate the phase interactions between oscillatory components and cross-frequency couplings. Some researchers (Vakorin, Mis, Krakovska, Bezgin, & McIntosh, 2013) compared the performance of the spectral, information-theoretic, and standard Granger causality in their relations to the observed phase differences and showed that the information-theoretic approach is the most effective against phase effects.

Relation of brain sites in some diseases such as schizophrenia and epilepsy has been studied. Granger causality (Cadotte, Demarse, Mareci, Parekh, Talathi, Huwang, Ditto, Ding, & Carney, 2010) was used to determine the direction of dynamic temporal relationships between LFP's recorded from bilateral microelectrode arrays before, during, and after seizures in an unprovoked spontaneous seizure model of temporal lobe epilepsy. The abnormal functional connectivity in schizophrenia by correlation between brain sites (Zalesky, Fornito, Egan, Pantelis, & Bullmore, 2011) was studied and altered coupling was demonstrated between regions and local decoherence within regions in the disease.

Temporal cross correlation of LFP's (Adhikari, Sigurdsson, Topiwala, & Gordon, 2010) is a simple effective method by means of which, the results were similar to the phase locking method (Sigurdsson, Stark, Karayiorgou, Gogos, & Gordon, 2010). In this study, we compared the cross correlation of potential activities (based onlyon LFP's) with phase locking (based on spikes and LFP's) through a simulation approach to reveal the advantages/disadvantages of the two methods from a processing neurophysiological point of view as a supplement to the previous studies.

2. Methods

2.1. Cross Correlation Method

Cross correlation of a couple of stochastic variables, X and Y, shown by r (X,Y) quantifies their linear coherence, i.e.

in which E[.] shows expectation function. Assuming X and Y as ergodic variables, it can be estimated by temporal realization of X and Y using the following equation:

in which N is the number of time samples available. As mentioned, this criterion shows linear coherence which means temporal relation during increase and decrease patterns.

2.2. Phase Locking Approach

Phase locking of a spike train to a local field potential signal is its coherence with the phase of the signal. As described by previous study (Sigurdsson, Stark, Karayiorgou, Gogos, & Gordon, 2010), considering the time when the kth spike occurs, a unit vector ‪(k) with a phase equal to the correspondent one at the signal is assigned, resulting in a set ofunit vectors (Nk is the number of spikes). If spikes are in agreement with some special phase of potential signal, they will gather near each other in complex space, while if they are not, they will fall around randomly. Defining sp as occurrence interval vector of multiunit spikes and phs as a vector containing phase of LFP during time, the phase lock power is evaluated by mean resultant length (MRL) of unit vectors assigned to spikes, i.e.

where ∣.∣ indicates absolute value.

2.3. Preprocessing

The potential signal was first passed from a band-pass filter to extract the information in the selected band (see next section on frequency band). The signal was then processed by Hilbert transform. The transform yields an imaginary signal which is the original one with a phase shift ofradians. A complex analytic signal was generated with the original signal as the real and its Hilbert transform as the imaginary part. The absolute value and the phase of the complex analytic signal in each time moment would, therefore, yield the instantaneous envelope and phase of the original band-pass filtered signal, respectively.

2.4. Estimation of Time Lag

Considering two brain sites, time lag of the first site in comparison with the second one is defined as the length of time after which it activates. This shows emission of processing activation through neural tissue.

2.4.1. Cross Correlation

Cross correlation is the simplest conceivable linear relation between two LFP signals. In order to estimate the time lag of LFP1 comparing LFP2, the cross correlation of LFP1 and LFP2 was taken for different time lags and the time lag respecting to the highest absolute value of cross correlation was then chosen as the result (Adhikari, Sigurdsson, Topiwala, & Gordon, 2010), illustrated by the following formula:

where ∣.∣ shows absolute value and LFP1 lag means LFP1 delayed by the amount of lag value. The absolute value indicates that the relation between the two signals may take opposite signs.

2.4.2. Phase Locking

In order to estimate the time lag of two brain sites by phase locking, power of the lock between the potential phases of one site, e.g. phase of LFP1 and multiunit spike train of the other was evaluated for a series of time lags. The result was the time lag by which the mentioned power would maximize (Adhikari, Sigurdsson, Topiwala, & Gordon, 2010), i.e.

in which phs1lag vector is the phase of LFP1 delayed by lag, and sp is a vector containing times at which spikes fired.

2.5. Simulation

Since the time lag between brain sites was not known on real brain data, we anticipated a simulation procedure by Matlab software to evaluate the capability of methods to detect the correct results. Signals of the first electrode in alpha band through all trials of real data were used as the first signal. The second signal was determined the same as the first one but with two optional predefined time lags of 0ms and 30ms comparing to it. To generate the spike train associated with the simulated second brain site, optionally, we put spikes on intervals when the phase of the first potential signal was in the range ofradians and shifted it with the amount of 0ms and 30ms. The methods of cross correlation and phase locking were expected to estimate the mentioned predefined time lag in the designed simulation.

2.6. Self Bias

In simulation, we did not shift the times spikes occurred comparing to the first signal but this was not certainly the fact in real data. There might be a time lag between spikes and potential of one brain site. This concept specially arises taking the fact into consideration that potentials and spikes, do not originate from common sources of activation. They are respectively attributed to dendritic and axonal bioelectric activities of neurons (Buzsáki, Anastassiou, & Koch, 2012). If there is a time lag between spikes and potential of one site, i.e. they are not synchronous; the lag should be neutralized by a time shift in opposite direction. This time shift is necessary since the time lag between a site comparing itself must be zero. Whether such self bias existed or not, could be answered using LFP and multiunit spike activity of the same brain sites.

3. Statistical Analysis

The available dataset was recorded from a gyrus on MT area of macaque brain which is a visual area that processes motion and velocity variations. A row of five electrodes was implanted on the gyrus. The electrodes were 0.3 mm apart from each other. At each electrode, LFP and spike signals were recorded simultaneously. Potential activity was recorded by a sampling frequency of 1000 Hz and spikes by 40 kHz. After signal was acquired, it was filtered, amplified and digitized. Omitting the beginning part before RDP (first three steps, Figure 1), each trial included 1000 samples equivalent to a 1-second data record.

Figure 1

Five steps in each trial data record, including fixation (by a plus sign), desired direction cue, fixation, RDP rep-resentation, random diversion in RDP direction.

The visual paradigm is illustrated in Figure 1. First by a plus sign, the monkey was fixated on the screen during the whole trial. Then a cue appeared that showed the location of the target to be attended. During each trial two moving random dot patterns (RDP) were shown on the screen, one inside and the other outside the receptive filed. Both patterns moved in one of eight possible directions selected randomly in each trial, i.e. 0 to 2radians withsteps (both had the same direction in each trial). At last, the monkey had to release a bar when the target RDP changed its direction.

4. Results

4.1. Frequency Band

Among frequency bands including Delta (0.5-4 Hz), Theta (4-8 Hz), Alpha (8-12 Hz), Beta (12-30Hz), Gamma1 (30-80Hz), and Gamma2 (80-150Hz), the one associated with maximum mutual information between the brain sites was considered (Figure 2). As a result, Alpha (8-12 Hz) was chosen which yielded meaningfully, the highest pairwise mutual information through electrodes comparing any other frequency band (wilcoxon ranksum test, p < 0.05). This means that the highest relation between the brain sites was seen in the mentioned band which was, therefore, the most proper one for lag estimation.

Figure 2

Mean mutual information (bits) of pair-electrode LFP's in different frequency bands with SEM. Alpha (8-12 Hz) demonstrated the highest mutual information in comparison with other band which was statistically meaningful (wilcoxon ranksum test, p < 0.05).

4.2. Simulation

We applied the simulation procedure on the local field potential signal recorded by the first electrode, containing129 trials. Estimation of time lags by (4) and (5) was done by a search in the interval of -50 ms to +50 ms with 5ms resolution. The results are shown in Table 1. Cross correlation and phase locking estimated the predefined synthetic time lags, i.e. -30 and 0, with no error. The results of phase locking were based on the linking 6-trial experiments, i.e. the spikes of every 6 trials (which contained 6 seconds of data record) together would yield a single result. This decision was made because using less-than-six combination of trials would cause error in lag estimation. As shown in Table 1, there was no error in phase locking method results with this consideration. Note that such combination was not seen necessary for cross correlation-based estimation of lags, since it yielded the predefined desired results by single 1 second-length trials.

Table 1

Estimated time lags (in ms) by cross correlation and phase locking for no delay and delayed simulation. The results of phase locking are based on 6-trial experiments. The results for both methods showed no errors throughout 129 trials.

	Correlation	Phase locking
Delay (-30ms)	-30	-30
No Delay	0	0

Frequency bands were investigated for minimum number of linked trials that were necessary to prevent the lag estimation error, depicted in Figure 3. As frequency increases, trials include a higher variety of fluctuations. In delta band, the errors would remain with any number of combined trials (even up to 129). In Theta, Alpha and Beta bands, 10, 6, and 2-trial experiments which included 10s, 6s, and 2s length of data were necessary to prevent errors. In Gamma1 and Gamma2 bands, there was no need to link trials which demonstrated there was enough variety of fluctuations in 1s trials in high frequency bands.

Figure 3

In each frequency band bars show minimum number of linked trials to form experiments with no error through 129 trials. In delta band estimation errors would remain by any number of linked trials (even up to 129). In Gamma1 and Gamma2 bands number of linked trials was 1, i.e. there was no need to link trials.

4.3. Self Bias

The simulation results implied that cross correlation would not contribute to self bias, i.e. a potential signal had no delay with itself. The same fact was implied for phase locking in the simulation results (see no delay results in Table 1).Table 2 shows the mean self bias of phase locking on real data for the five electrode records through 129 trials with a 2ms resolution, i.e. how long after potential activity, spikes would fire. In all cases, the lag was positive which means that spikes fired later than local field potential activity. Second and third electrodes suffered from self biases meaningfully greater than zero but it was not statistically proved for other electrodes (signed rank wilcoxon test, p < 0.05). Since spike and LFP are respectively attributed to axonal and dendritic bioelectric activity (Buzsáki, Anastassiou, & Koch, 2012), from a neurophysiological point of view, time latency of spike train interpreted as tissue output in regard to LFP interpreted as tissue input is expected. Hence, findings reported in Table 2 are in agreement with the suggested origin of spike and LFP signals. As a solution to the problem, the bias must be neutralized by a time shift of spike occurrence moments in opposite direction before phase locking being used for lag estimation between the brain sites. This time shift is necessary since the lag between a site comparing itself must be zero.

Table 2

Self bias of five electrodes and significance of rejection of the null hypothesis that the self bias had zero median (signed rank wilcoxon test, p < 0.05).

Electrode	1	2	3	4	5
Self Bias (ms)	4.0	6.7	5.5	3.5	2.5
Significance	0.066	p < 0.05	p < 0.05	0.155	0.320

5. Discussion

We compared the capability of two methods, i.e. cross correlation (based on LFP) and phase locking (based on LFP and spikes) to estimate the time lags of brain sites on MT area of macaque brain and as a result, the former demonstrated higher efficiency. Cross correlation estimated the true time lag with no error through129 trials, but phase locking estimated the true results with no error by 6-trial linked experiments. This showed, for the first time, that phase locking estimation of time lag especially in lower frequency bands was highly dependent on data length and therefore, not time efficient. Hence, trials must contain at least a minimum variety of signal fluctuations to get valid results from the method. The threshold would increase in lower frequency bands which meant enough variety of fluctuations entailed longer data length in low frequency bands. In addition, we introduced for the first time, the problem of self bias of phase locking method estimating the time lag between brain sites, which was neurophysiologically in agreement with the neuronal origins of spike and LFP signals. The bias was estimatedms through 5 recording sites in our data set. As a solution to the self bias problem, it must be removed by a time shift on spike occurrence moments, otherwise, non-zero time lag between a site with itself will be reported which is obviously, a conflict. Since the bias is not the same for different sites, it should be done as a preprocessing step on each recording site, so that the spike train and local field potential become synchronous. In previous studies (Sigurdsson, Stark, Karayiorgou, Gogos, & Gordon, 2010; Adhikari, Sigurdsson, Topiwala, & Gordon, 2010), this problem led to biases in the results. Applying this correction e.g. on the latter study, the negative time lag between the two brain regions (-24.5ms) will be summed by some positive scalar, which leads to a decrease in absolute lag value reported between two brain regions.

Importantly, as no-delay results of simulation showed, the self bias was not caused by phase locking method itself. The bias resulted from the asynchrony of spikes and LFP's of the same site, hence, it had neurophysiological basis, not a signal processing one.