elephant.causality.granger.pairwise_granger¶
- elephant.causality.granger.pairwise_granger(signals, max_order, information_criterion='aic')[source]¶
Determine Granger Causality of two time series
- Parameters:
- signals(N, 2) np.ndarray or neo.AnalogSignal
A matrix with two time series (second dimension) that have N time points (first dimension).
- max_orderint
Maximal order of autoregressive model.
- information_criterion{‘aic’, ‘bic’}, optional
A function to compute the information criterion:
bic for Bayesian information_criterion,
aic for Akaike information criterion,
Default: ‘aic’
- Returns:
- Causality
- A namedtuple with the following attributes:
- directional_causality_x_yfloat
The Granger causality value for X influence onto Y.
- directional_causality_y_xfloat
The Granger causality value for Y influence onto X.
- instantaneous_causalityfloat
The remaining channel interdependence not accounted for by the directional causalities (e.g. shared input to X and Y).
- total_interdependencefloat
The sum of the former three metrics. It measures the dependence of X and Y. If the total interdependence is positive, X and Y are not independent.
- Denote covariance matrix of signals
X by C|X - a real number Y by C|Y - a real number (X,Y) by C|XY - a (2 times 2) matrix
- directional causality X -> Y given by
log(C|X / C|XY_00)
- directional causality Y -> X given by
log(C|Y / C|XY_11)
- instantaneous causality of X,Y given by
log(C|XY_00 / C|XY_11)
- total interdependence of X,Y given by
log( {C|X cdot C|Y} / det{C|XY} )
- Raises:
- ValueError
If the provided signal does not have a shape of Nx2.
If the determinant of the prediction error covariance matrix is not positive.
- Warns:
- UserWarning
If the log determinant of the prediction error covariance matrix is below the tolerance level of 1e-7.
Notes
The formulas used in this implementation follows (Ding et al., 2006). The only difference being that we change the equation 47 in the following way: -R(k) - A(1)R(k - 1) - … - A(m)R(k - m) = 0. This forumlation allows for the usage of R values without transposition (i.e. directly) in equation 48.
Examples
Example 1. Independent variables.
>>> import numpy as np >>> from elephant.causality.granger import pairwise_granger >>> pairwise_granger(np.random.uniform(size=(1000, 2)), max_order=2) # noqa Causality(directional_causality_x_y=-0.0, directional_causality_y_x=0.0, instantaneous_causality=0.0, total_interdependence=0.0)
Example 2. Dependent variables. Y depends on X but not vice versa.
\[\begin{split}\begin{array}{ll} X_t \sim \mathcal{N}(0, 1) \\ Y_t = 3.5 \cdot X_{t-1} + \epsilon, \; \epsilon \sim\mathcal{N}(0, 1) \end{array}\end{split}\]In this case, the directional causality is non-zero.
>>> np.random.seed(31) >>> x = np.random.randn(1001) >>> y = 3.5 * x[:-1] + np.random.randn(1000) >>> signals = np.array([x[1:], y]).T # N x 2 matrix >>> pairwise_granger(signals, max_order=1) # noqa Causality(directional_causality_x_y=2.64, directional_causality_y_x=-0.0, instantaneous_causality=0.0, total_interdependence=2.64)