pyhctsa.operations.correlation.joint_non_gaussianity

pyhctsa.operations.correlation.joint_non_gaussianity(y, tau='ac', m=2, theiler_win=('ac', 1), max_n=10000)

Tests for non-Gaussianity of the joint, time-lagged embedding distribution.

Port of hctsa’s CO_JointNonGaussianity. Embeds the time series in m dimensions at time delay tau (e.g., the pair (x_t, x_{t+tau}) for m=2, or the triple (x_t, x_{t+tau}, x_{t+2tau}) for m=3) and tests whether the resulting point cloud is consistent with a multivariate Gaussian.

A linear (e.g., AR(1)) Gaussian process has a Gaussian marginal and a Gaussian joint embedding distribution; a nonlinear or non-reversible process can look Gaussian marginally while its lagged joint distribution is visibly non-elliptical (curved, multimodal, or heavy/light-tailed along directions the marginal alone cannot see). This tests the whole joint shape rather than one moment combination (cf. trev() and tc3()).

Two complementary statistics are based on Mardia’s classical multivariate normality measures [1], chosen because they generalize to any embedding dimension via the same formula and reduce, at m=1, to ordinary skewness/kurtosis:

  • Mardia’s multivariate skewness, b1: detects asymmetry/curvature of the joint distribution (e.g., a banana-shaped point cloud). Its population value is 0 for any joint Gaussian.

  • Mardia’s multivariate kurtosis, b2: detects joint tail weight/peakedness relative to a Gaussian ellipsoid. Its population value is m(m+2) for any joint Gaussian (8 at m=2, 15 at m=3).

As a complementary check, the squared Mahalanobis distances of each embedded point to the sample mean (the per-point terms underlying Mardia’s kurtosis) are compared with their theoretical distribution under joint Gaussianity, chi^2_m, via a Kolmogorov-Smirnov D-statistic. This can catch departures (e.g., a bimodal or ring-shaped cloud) that the two summary moments can miss.

The skewness statistic excludes near-diagonal pairs (|i-j| <= theiler_win) from its double sum: for correlated jointly-Gaussian points, the third moment of their Mahalanobis inner product is not zero, so nearby, strongly autocorrelated pairs bias the raw statistic away from zero even under true joint Gaussianity. This removes most, but not all, of the bias (the residual is the classical small-sample bias of whitening with the sample covariance of the points being tested).

References

Parameters:
y : array-like

The input time series.

tau : int or {'ac', 'ac1e', 'mi'}, optional

The time delay for the embedding (see get_tau()). Default 'ac'.

m : int, optional

The embedding dimension (default 2, the pairwise joint distribution (x_t, x_{t+tau}); set to 3 for the triple-wise joint distribution). hctsa’s false-nearest-neighbors choice of m ({'fnn', th}) is not available.

theiler_win : (str, float) or int or float, optional

The number of temporally adjacent embedded points excluded from the skewness double sum (|i-j| <= theiler_win): ('ac', k) for k times the first zero-crossing of the autocorrelation function, ('ac1e', k) for k times its 1/e delay, or a number of samples (see theiler_window()). Default ('ac', 1).

max_n : int or 'full', optional

The maximum number of embedded points used for the skewness statistic (default 10000; 'full' to disable). Only the first max_n embedded points enter the skewness, as in hctsa (where the cap bounds an N x N Gram matrix). Here the skewness double sum is evaluated exactly through the third-moment tensor in O(N d^3), so the cap does not save time; it keeps values consistent with hctsa. The mean, covariance, kurtosis and KS statistic always use the full embedded series.

Returns:

  • mardiaSkew: Mardia’s raw multivariate skewness (Theiler-windowed);

  • mardiaKurt: Mardia’s raw multivariate kurtosis;

  • mahalKSstat: the Kolmogorov-Smirnov D-statistic of the squared Mahalanobis distances against chi^2_m.

All are unitless departure-from-joint-Gaussianity magnitudes with no attached significance level (see Notes). The output is a single NaN if the embedding fails, has too few points, or has a near-singular covariance (mardiaSkew alone is NaN if the Theiler window leaves no pairs).

Return type:

dict or float

Notes

Only raw statistics are returned, not p-values. Mardia’s classical asymptotic null distributions assume the N embedded points are iid draws, but consecutive embedded vectors overlap in m-1 coordinates and are strongly autocorrelated, which inflates the naive test statistics (giving far more false positives than the nominal level even for a linear-Gaussian AR(1) process, and more so at higher m). For significance testing, compare to the distribution over surrogates.

The KS statistic is computed against the exact chi^2_m CDF; hctsa evaluates it on a table of the CDF at values rounded to 1e-6, so the two agree to about 1e-6.