elephant.signal_processing.cross_correlation_function

elephant.signal_processing.cross_correlation_function(signal, channel_pairs, hilbert_envelope=False, n_lags=None, scaleopt='unbiased')

Computes an estimator of the cross-correlation function (Stoica et al., 2005).

\[\begin{split}R(\tau) = \frac{1}{N-|k|} R'(\tau) \\\end{split}\]

where \(R'(\tau) = \left<x(t)y(t+\tau)\right>\) in a pairwise manner, i.e.:

signal[channel_pairs[0,0]] vs signal[channel_pairs[0,1]],

signal[channel_pairs[1,0]] vs signal[channel_pairs[1,1]],

and so on.

The input time series are z-scored beforehand. scaleopt controls the choice of \(R_{xy}(\tau)\) normalizer. Alternatively, returns the Hilbert envelope of \(R_{xy}(\tau)\), which is useful to determine the correlation length of oscillatory signals.

Parameters:

signal(nt, nch) neo.AnalogSignal

Signal with nt number of samples that contains nch LFP channels.

channel_pairslist or (n, 2) np.ndarray

List with n channel pairs for which to compute cross-correlation. Each element of the list must contain 2 channel indices. If np.ndarray, the second axis must have dimension 2.

hilbert_envelopebool, optional

If True, returns the Hilbert envelope of cross-correlation function result. Default: False

n_lagsint, optional

Defines the number of lags for cross-correlation function. If a float is passed, it will be rounded to the nearest integer. Number of samples of output is 2*n_lags+1. If None, the number of samples of the output is equal to the number of samples of the input signal (namely nt). Default: None

scaleopt{‘none’, ‘biased’, ‘unbiased’, ‘normalized’, ‘coeff’}, optional

Normalization option, equivalent to matlab xcorr(…, scaleopt). Specified as one of the following.

• ‘none’: raw, unscaled cross-correlation

\[R_{xy}(\tau)\]

• ‘biased’: biased estimate of the cross-correlation:

\[R_{xy,biased}(\tau) = \frac{1}{N} R_{xy}(\tau)\]

• ‘unbiased’: unbiased estimate of the cross-correlation:

\[R_{xy,unbiased}(\tau) = \frac{1}{N-\tau} R_{xy}(\tau)\]

• ‘normalized’ or ‘coeff’: normalizes the sequence so that the autocorrelations at zero lag equal 1:

\[R_{xy,coeff}(\tau) = \frac{1}{\sqrt{R_{xx}(0) R_{yy}(0)}} R_{xy}(\tau)\]

Default: ‘unbiased’

Returns:

cross_corrneo.AnalogSignal

Shape: [2*n_lags+1, n] Pairwise cross-correlation functions for channel pairs given by channel_pairs. If hilbert_envelope is True, the output is the Hilbert envelope of the pairwise cross-correlation function. This is helpful to compute the correlation length for oscillating cross-correlation functions.

Raises:

ValueError

If input signal is not a neo.AnalogSignal.

If channel_pairs is not a list of channel pair indices with shape (n,2).

If hilbert_envelope is not a boolean.

If n_lags is not a positive integer.

If scaleopt is not one of the predefined above keywords.

Examples

import neo
import numpy as np
import quantities as pq
import matplotlib.pyplot as plt
from elephant.signal_processing import cross_correlation_function
dt = 0.02
N = 2018
f = 0.5
t = np.arange(N)*dt
x = np.zeros((N,2))
x[:,0] = 0.2 * np.sin(2.*np.pi*f*t)
x[:,1] = 5.3 * np.cos(2.*np.pi*f*t)

# Generate neo.AnalogSignals from x and find cross-correlation

signal = neo.AnalogSignal(x, units='mV', t_start=0.*pq.ms,
         sampling_rate=1/dt*pq.Hz, dtype=float)
rho = cross_correlation_function(signal, [0,1], n_lags=150)
env = cross_correlation_function(signal, [0,1], n_lags=150,
    hilbert_envelope=True)

plt.plot(rho.times, rho)
plt.plot(env.times, env) # should be equal to one
plt.show()

(Source code, png, hires.png, pdf)