Probing high order dependencies with information theory

TL;DR

This paper introduces the use of Auto Mutual Information and entropy rate to detect high-order dependencies in signals, revealing details missed by traditional correlation analysis and distinguishing between complex non-Gaussian processes.

Contribution

It demonstrates how Auto Mutual Information and entropy rate can assess refined, higher-order temporal dependencies in signals, surpassing traditional correlation-based methods.

Findings

Auto Mutual Information captures higher-order dependencies.

Entropy rate distinguishes non-Gaussian processes with similar covariance.

Estimators' performance analyzed via Monte Carlo simulations.

Abstract

Information theoretic measures (entropies, entropy rates, mutual information) are nowadays commonly used in statistical signal processing for real-world data analysis. The present work proposes the use of Auto Mutual Information (Mutual Information between subsets of the same signal) and entropy rate as powerful tools to assess refined dependencies of any order in signal temporal dynamics. Notably, it is shown how two-point Auto Mutual Information and entropy rate unveil information conveyed by higher order statistic and thus capture details of temporal dynamics that are overlooked by the (two-point) correlation function. Statistical performance of relevant estimators for Auto Mutual Information and entropy rate are studied numerically, by means of Monte Carlo simulations, as functions of sample size, dependence structures and hyper parameters that enter their definition. Further, it is shown how Auto Mutual Information permits to discriminate between several different non Gaussian processes, having exactly the same marginal distribution and covariance function. Assessing higher order statistics via multipoint Auto Mutual Information is also shown to unveil the global dependence structure fo these processes, indicating that one of the non Gaussian actually has temporal dynamics that ressembles that of a Gaussian process with same covariance while the other does not.

Tables

Table 1. TABLE I: Statistical test. p 𝑝 p -value of Wilcoxon test between the statistics of I ^ ( 1 , 3 ) ​ ( 1 ) superscript ^ 𝐼 1 3 1 \hat{I}^{(1,3)}(1) and I ^ ( 2 , 2 ) ​ ( 1 ) superscript ^ 𝐼 2 2 1 \hat{I}^{(2,2)}(1) .

	$Y_{𝕀}$	$Y_{F_1}$	$Y_{F_2}$
Short range correlations	0.34	0.96	0.72
Long range correlations	$6e^{-15}$	$5e^{-20}$	$2e^{-3}$

Equations

$${\text{x}}_t^{(m,\tau )}=\left(x_t,x_{t-\tau },\cdots ,x_{t-(m-1)\tau }\right),$$

$$H({\text{x}}_t^{(m,\tau )})$$

$$h^{(m,\tau )}(x_t)$$

$$=H^{(m+1,\tau )}(X)-H^{(m,\tau )}(X)$$

$$=h^{(m,\tau )}(X).$$

$$I({\text{x}}_t^{(m,\tau )},$$

$$\ln \left(\frac{p\left({\text{x}}_t^{(m,\tau )},{\text{y}}_{t^{\prime }}^{(p,\tau )}\right)}{p({\text{x}}_t^{(m,\tau )})p({\text{y}}_{t^{\prime }}^{(p,\tau )})}\right)\mathrm{d}{\text{y}}_{t^{\prime }}^{(p,\tau )}\mathrm{d}{\text{x}}_t^{(m,\tau )}$$

$$I^{(m,p,\tau )}(X)=I({\text{x}}_t^{(m,\tau )},{\text{x}}_{t+p\tau }^{(p,\tau )}).$$

$$\left({\text{x}}_t^{(m,\tau )},{\text{x}}_{t+p\tau }^{(p,\tau )}\right)={\text{x}}_{t+p\tau }^{(m+p,\tau )},$$

$$h^{(m,\tau )}(X)=H(X)-I^{(m,1,\tau )}(X),$$

$$H(X)=$$

$$I^{(m,p,\tau )}(X)=$$

$$h^{(m,\tau )}(X)=$$

$$\hat{H}(X)=-\psi (k)+\psi (N)+\frac{d}{N}\sum _i^N\log ({ϵ}_i(k))$$

$$\hat{I}(X,Y)=\psi (k)+\psi (N)-⟨\psi (n_x)+\psi (n_y)⟩$$

$$H(Y)$$

$$=\frac{1}{2}\ln (2\pi e{\sigma }_X^2)+⟨\ln F{⟩}_X$$

$$I^{(m,p,\tau )}(Y)$$

$$=\frac{1}{2}\ln \left(\frac{|{\Sigma }_X^{(m)}||{\Sigma }_X^{(p)}|}{|{\Sigma }_X^{(m+p)}|}\right)$$

$$h^{(m,\tau )}(Y)=$$

$$F_1(x)=F_Y^{-1}\circ \Phi (x),$$

$$F_2(x)=F_Y^{-1}\circ \left(2(\Phi (|x|)-1/2)\right),$$

$${\mu }_{X_1}$$

$$c_{X_1}(\tau )=\frac{1}{{\sigma }_{X_1}^2}\ln ((e^{{\sigma }_{X_1}^2}-1)c_Y(\tau )+1).$$

$$c_Y(\tau )$$

Full text

I Introduction

Context. To characterize a (stationary) process $X(t)$, evolving along a given dimension $t$, e.g., time, the precise assessment of both its static properties (marginal distribution or one-point statistics) and temporal dynamics (joint distributions or multi-point statistics) is crucial.

Correlation and spectral estimates, providing practitioners with valuable insight into temporal dynamics are key tools in statistics and data analysis and have been widely used in many different fields. In neurosciences, the analysis of correlation between neuron firings [1, 2] or between stimulus and brain activity [3, 4] is essential. In fluid mechanics, Kolmogorov’s empirical theory of turbulence [5] characterizes the multiscale distribution of energy with the power spectrum of the velocity field [6]. The analysis of correlations also supports useful descriptions in terms of complex networks for, e.g., ecology [7], climatology [8] or traffic [9].

It is however well known that, while the correlation function, or equivalently the power spectral density, exhaust the description of joint statistics and hence of temporal dynamics for jointly Gaussian processes, they provide characterizations of temporal dynamics for general processes, that are limited in at least two respects: not only they are restricted to two-point joint statistics but also to a specific measure of these two point-statistics. In general, the full analysis of temporal dynamics require the use of multi-point, or higher order, joint statistics. The goal of the present work is to address such an issue using information theoretic measures.

Related work. To assess dependencies beyond two-point correlations, multipoint correlation functions were defined. For instance, three-point correlation functions have been used in cosmology, where they were related to cosmological parameters [10, 11]. Alternatively, higher order-spectra, that combine Fourier transforms at several different frequencies and therefore takes into account higher order statistics, have also been used to generalize the power spectrum beyond second-order statistics [12, 13]. For instance, bi-spectrum and tri-spectrum were analyzed in depth in e.g., [14, 15] and used in applications as different as cosmology [16, 17], or medicine [18, 19] and other fields [20, 21]. While of popular use in applications, higher-order correlation functions and spectra suffers from major caveats, mostly related to estimation issues. Indeed, the assessment of order $m$-statistics with higher-order correlation functions and spectra requires the estimation of $m$-time point statistics, which, for a given fixed sample size, can only yield decreasing estimation performance when order $m$ is increased.

Alternatively, Shannon’s information theory, and the corresponding information theoretic measures, entropy, entropy rate and mutual information, can be considered as an efficient way to characterize high order statistics [22, 23, 24]. Applications of entropy rate were yet reported in a wide range of domains, such as biology [25], dynamical systems [26], fluid turbulence [27] or analysis of languages [28]. Mutual information is also very popular, especially in computer sciences, e.g., in machine learning [29], and other very different domains [30, 31, 32, 33, 34, 35]. Auto-Mutual Information (AMI) — sometimes called Information Storage [36] —, defined as the mutual information between the analyzed process and a delayed version of itself, has been used in different domains [37, 38, 39, 40]. Information theory additionally provides tools to explore causality relationships between processes [41, 42]. However, within this paper we are mainly interested in the characterization of their dependence structures.

Entropy rate and AMI have been envisaged as tools to describe temporal dynamics in time series beyond second-order statistics, or even at second-order statistics, in a richer way than correlation or Fourier spectrum, see  [43, 37, 39, 27].

Goals, contributions and outline. The main goals of the article are i) to show the ability of AMI and entropy rate to describe high order dependences, even when considering only two-point interactions. Even though this property of information theory tools is well known, we are able to show it by analyzing processes with identical correlation and marginal distribution but different high order dependences. AMI and entropy rate appear as straightforward generalizations of correlation function and power spectrum when considering non-Gaussian processes analysis. ii) to show the effect of analysing interactions between more than two points, combining Takens embedding [44] with AMI. We characterize the complexity of the dependence structure of a process by measuring AMI for embedding dimension higher than 2. More precisely, we analyze the differences between power law decay and exponential decay dependence structures.

To demonstrate the ability of information theory to describe high order statistics, we illustrate the use of entropy rate and auto-mutual information to characterize non Gaussian processes across scales. As an example, we consider two synthetic log-normal processes with different dependence structures but with identical correlation function and identical marginal statistics. In section II, we introduce and define the theoretical information theory measures and our estimators. In section III, we interpret the information theory tools used along the article, we present the processes with identical correlation function and one-point statistics but different dependence structure that we study and finally we illustrate how AMI and entropy rate allow to distinguish between them. In section IV-A, we report the convergence and the performance of our estimators of AMI and entropy rate and their ability to measure high order dependencies. We present their evolution across scales and emphasized that, contrary to classical correlation analysis, they allow a fine characterization of non Gaussian processes.

II Entropy, entropy rate and auto mutual information

II-A Definitions

To analyze the temporal dynamics of univariate stationary processes $X=\left\{x_{t}\right\}_{t\in\mathbb{R}}$, it has been proposed in [44] to rely on the time delay-embedding procedure. It amounts to constructing the $m$-dimensional process $X^{(m,\tau)}$, whose elements are $m$-dimensional vectors:

[TABLE]

where the delay $\tau$ implicitly defines an analysis scale and where the embbeding dimension $m$ controls the order of the statistics of $X$ which are analyzed. The information theoretic measures defined below aim to characterize temporal dynamics via the analysis of the $m$-time point joint statistical distributions $p(\textbf{x}^{(m,\tau)}_{t})$.

Within this article, we study ergodic stationary processes and then $p(\textbf{x}^{(m,\tau)}_{t})$ depend on $m$ and $\tau$ but not on $t$. Consequently $p(\textbf{x}^{(m,\tau)}_{t})=p(X^{(m,\tau)})$.

Shannon entropy, is defined as a functional of $p(\textbf{x}^{(m,\tau)}_{t})$ [22]:

[TABLE]

By nature, it depends on the full joint distribution and therefore on statistics of any order. However, for stationary processes it does not depend on the mean value of the marginal distribution. For $m=1$, Shannon entropy does not relate to temporal dynamics but to one-point statistics and thus obviously does not depend on $\tau$ and is hence denoted $H(X)$.

Entropy rate of order $m$, $h^{(m,\tau)}$, quantify variations of Shannon entropy between two successive time-embedded versions of a stationary process [22, 24]:

[TABLE]

Entropy rate $h^{(m,\tau)}$, like $H^{(m,\tau)}$, probes joint statistics of any order, but in addition quantifies how much information is gained in the analysis of temporal dynamics by increasing by one the order of the embedding.

Mutual information $I(\textbf{x}_{t}^{(m,\tau)},\textbf{y}_{t^{\prime}}^{(p,\tau)})$ quantifies the amount of information shared by two processes $X$ and $Y$, via their embbeded versions $\textbf{x}_{t}^{(m,\tau)}$ and $\textbf{y}_{t^{\prime}}^{(p,\tau)}$ [22, 24]:

[TABLE]

When $X$ and $Y$ are jointly stationary processes, $I(\textbf{x}_{t}^{(m,\tau)},\textbf{y}_{t^{\prime}}^{(p,\tau)})=I(X_{t}^{(m,\tau)},Y_{t^{\prime}}^{(p,\tau)})$ obviously depends on $t^{\prime}-t$ only.

Mutual information can be read as the Kullback-Leibler divergence between the joint probability $p(\textbf{x}_{t}^{(m,\tau)},\textbf{y}_{t^{\prime}}^{(p,\tau)})$ and its counterpart under independence hypothesis $p(\textbf{x}_{t}^{(m,\tau)})p(\textbf{y}_{t^{\prime}}^{(p,\tau)})$.

To analyze the higher statistical order temporal dynamics (also referred to as nonlinear temporal dynamics) of a single stationary process $X$, one can use auto-mutual information [24, 36], defined as the mutual information between $\textbf{x}_{t}^{(m,\tau)}$ and $\textbf{x}_{t+p\tau}^{(p,\tau)}$:

[TABLE]

Because the concatenation of $\textbf{x}_{t+p\tau}^{(p,\tau)}$ and $\textbf{x}_{t}^{(m,\tau)}$ corresponds exactly to the $(m+p)$-dimensional time-embedded vector $\textbf{x}_{t+p\tau}^{(m+p)}$:

[TABLE]

$I^{(m,p,\tau)}(X)$ quantifies the shared information between the past $m$-point dynamics and the future $p$-point dynamics. Auto mutual information quantifies the deviation from independence between the two arguments, beyond the simpler decorrelation.

These information theoretic quantities all depend on the chosen time lag $\tau$. When $\tau\to 0$, dependences between two consecutive coordinates of the embedding vector increase, therefore inducing an increase in $I^{(m,1,\tau)}(X)$ and a decrease in $h^{(m,\tau)}(X)$. Conversely, when $\tau\to+\infty$, dependencies decrease to [math], therefore inducing an increase of $h^{(m,\tau)}(X)$ up to the Shannon entropy $H(X)$ of the process, its largest value, and a decrease of $I^{(m,1,\tau)}(X)\rightarrow 0$.

Furthermore, entropy rate can then be expressed as [27]:

[TABLE]

which shows that there are two different contributions to the entropy rate: First, $H(X)$, which only depends on the one-point distribution and hence is a static property ; second, $I^{(m,1,\tau)}(X)$ which gathers all information conveyed by linear and non linear temporal dynamics, irrespective of the variance of the distribution [27].

II-B Gaussian process

When $X$ is a stationary jointly Gaussian process, hence fully defined by its variance $\sigma_{x}^{2}$ and normalized correlation function $c_{X}(\tau)$, analytical expressions for entropy, entropy rate and auto mutual information are available [45]:

[TABLE]

with $\Sigma^{(m)}$ the $m$-dimensional correlation matrix, and $||$ denoting the determinant. By definition, $\Sigma^{(1)}=1$.

II-C Estimation procedures

To estimate these information theoretic quantities, two main categories of non-parametric procedures are available, based on either kernel density estimations (KDE) or nearest neighbor search (NNS). While the former essentially counts the number $k(\epsilon)$ of the observed samples in interval of a priori set size $\epsilon$, the latter rather measure the size $\epsilon(k)$ of box that contain an a priori set number of $k$ samples.

Because this second strategy has been documented as providing low bias estimates [46, 47], we chose $k$-nearest neighbors procedures to estimate the entropy $H(X)$ [48, 49, 50]:

[TABLE]

with $\psi$ the Digamma function, $d$ the dimension of $X$, $N$ the number of samples of $X$, and using the Chebyshev distance ($L_{\infty}$ metric).

To compute the mutual information, we use the $k$-nearest neighbors procedure presented by Kraskov, et al. [46]:

[TABLE]

The procedure chooses $k$ for the joint distribution $Z(X,Y)$ and hence obtains the distance $\epsilon$ to the $k$th nearest neighbour. Then using $\epsilon$, it calculates the number of neighbours $n_{x}+1$ lying on a segment $\textbf{x}=x\pm\frac{\epsilon}{2}$ of the X dimension, i.e. lying on $-\infty\leq Y\leq\infty$, $x-\epsilon/2\leq X\leq x+\epsilon/2$. Consequently, $\epsilon/2$ is the distance to the $n_{x}+1$st neighbour of $x$. The procedure uses the same development for $Y$ and obtains that $\epsilon/2$ is the distance to the $n_{y}+1$st neighbour of $y$.

Because benefits against a direct use of eq.(2) were documented in [51], $h^{(m,\tau)}(X)$ is estimated plugging the above estimates into eq.(5).

For an stationary process $X$, $p(\textbf{x}^{(m,\tau)}_{t})$ does not depend on $t$, and consequently, we can exploit this property to estimate the information theoretic quantities of $X$ by averaging on all the temporal samples of the process.

To perform an analysis across scales by varying the parameter $\tau$ we use an adapted Theiler prescription [52] to eliminate the contribution of dependencies from scales smaller than $\tau$.

III Temporal dynamics and dependence structure

The aim of this section is to try to characterize what part of the statistics of a given stationary process $Y$ contribute to the values taken by information theoretic measures $H$, $h^{(m,\tau)}$ and $I^{(m,p,\tau)}$. To that end, we conduct a formal study of the specific yet very instructive case where processes $Y$ with correlation function $c_{Y}(\tau)$ and marginal distributions $p_{Y}$ are obtained as point bijective transformation $F$, $Y=F(X)$ of an underlying Gaussian process $X$ with correlation function $c_{X}$, designed to reach the targeted covariance $c_{Y}$. Processes $Y$ constructed in this way depend on $F$ and $c_{X}(\tau)$ [53, 54].

III-A Analytical calculations and interpretations

Given any bijective transformation $F$ that maps a Gaussian signal $X$ into a signal $Y=F(X)$, the information theoretic quantities of $Y$ can be computed as:

[TABLE]

with $\langle\ln F\rangle_{X}=\int p(\textbf{x}_{t}^{(1)})\ln F(x)dx$, hence a constant, $\sigma_{X}^{2}$ and $\Sigma_{X}^{(m)}$ are the variance and the $m$-dimensional covariance matrix of the underlying Gaussian process $X$.

These analytical calculations permit to yield understanding of general validity regarding the nature of the information conveyed by each of these three information theoretic quantities:

i) Entropy $H$ is only related to the marginal distribution of the process $p_{Y}$, hence is quantifying static properties only.

ii) Entropy rate $h^{(m,\tau)}$ mixes up marginal distribution features (variance and pointwise transform $F$), hence static properties with dynamical properties and dependencies of any order, as indeed $\Sigma_{X}^{(m)}$ depends jointly on $c_{Y}$ and on $F$, hence on the joint distribution of $\textbf{y}_{t}^{(m,\tau)}$. The fact that $h^{(m,\tau)}$ gather both static and dynamics aspects of time series likely explains why it has been observed in many applications to be a often discriminant feature, compared to other features [40]. However, the entanglement of static and dynamic properties may also be considered a drawback for the analysis, characterization and understanding of data properties.

iii) Interestingly, auto mutual information $I^{(m,p,\tau)}$ performs a characterization mainly related to the dependence structure of $Y$, as indeed $\Sigma_{X}^{(m)}$ depends on the joint distribution of $\textbf{y}_{t}^{(m,\tau)}$. It therefore quantifies temporal dynamics, irrespective of static properties. Thus, AMI is the sole measure that focuses exclusively on higher order dynamics and on the joint dependence structure, beyond correlation.

Often in the scientific literature, $h^{(m,\tau)}$ and $I^{(m,p,\tau)}$ are referred to as nonlinear features, to emphasize that they characterize temporal dynamics beyond the mere use of the correlation functions or power spectral densities. This is indeed correct. Even when two-point statistics are analyzed, i.e., $m=p=1$, $I^{(1,1,\tau)}$ reflects dependencies, including and encompassing correlations as the whole joint two-point distribution function is involved. Nevertheless, $I^{(m,p,\tau)}$ remains mostly driven by the correlation function $c_{Y}$, and hence by linear effects. This is also the case for $h^{(m,\tau)}$, driven mostly both by the static variance $\sigma^{2}_{Y}$ and dynamical correlation function $c_{Y}$.

III-B Process definition

III-B1 Analytical derivations

To further illustrate that indeed $I^{(m,p,\tau)}$ and $h^{(m,\tau)}$ may capture differences in any-order dependencies beyond correlations, let us study the cases where two stationary processes $Y_{F_{1}}$ and $Y_{F_{2}}$ have exact same marginal distribution $p_{Y}$ and exact same covariance $c_{Y}$, but are yet obtained from two different point transformations $F_{1}$ and $F_{2}$ — applied to Gaussian processes $X_{1}$ and $X_{2}$ — and hence have different joint distributions and hence different dependence structures.

Let denote $F_{Y}$ the cumulative distribution function associated to the marginal distribution $p_{Y}$. A natural candidate to obtain a process with marginal $p_{Y}$ is to apply the point transform

[TABLE]

to a Gaussian process $X_{1}$, where $F_{Y}$ is the cumulative distribution function of the marginal distribution $p_{Y}$ and $\Phi(x)$ the cumulative distribution function of a centered standard Gaussian variable, and where $X_{1}$ is a jointly Gaussian process with covariance function $c_{X_{1}}(\tau)$ to reach the targeted $c_{Y}(\tau)$ [53, 54].

However, many other transforms can yield processes $Y$ with the exact same marginal distribution $p_{Y}$ and the exact same covariance $c_{Y}$. In [54], it was for example proposed to use

[TABLE]

on a jointly Gaussian process $X_{2}$. This thus yields a radically different joint dependence structure for $Y_{F_{2}}$, despite having the same prescribed correlation function as $Y_{F_{1}}$. As a consequence, $Y_{F_{1}}$ and $Y_{F_{2}}$ have different AMI functions $I^{(m,p,\tau)}(Y_{F_{1}})\neq I^{(m,p,\tau)}(Y_{F_{2}})$, and entropy rates, and hence different joint dependence structure and thus temporal dynamics. Conversely, their entropies are identical $H(Y_{F_{1}})=H(Y_{F_{2}})$. This is not surprising , as the entropy is completely defined by the one-point distribution.

III-B2 The log-normal example

To quantify the qualitative analysis proposed in Section III-A , let us consider stationary processes $Y_{F_{1}}$ and $Y_{F_{2}}$ having the same log-normal marginal distribution $p_{Y}$, chosen as an example of distributions that may complicate estimation because of its being a frontier between distributions with all moments finite and heavy-tail distributions.

For process $Y_{F_{1}}$ obtained with the natural transform (14), $F_{1}$ simplifies to $F_{1}(x)=\exp(x)$, which is a bijective transformation. One can therefore use eqs.(11) to (13) above with $\langle\ln F_{1}\rangle_{X_{1}}=\mu_{X_{1}}$, $\sigma_{X_{1}}$ and $\Sigma_{X_{1}}^{(m)}$ being the mean, variance and $m$-dimensional covariance matrix of the Gaussian variable $X_{1}$ [45]. The entropy $H(Y_{F_{1}})$ can be computed using the mean $\mu_{Y}$ and variance $\sigma_{Y}^{2}$ of $Y_{F_{1}}$:

[TABLE]

Writing $I^{(m,p,\tau)}(Y_{F_{1}})$ thus only requires to express the correlation function of the underlying Gaussian process $X_{1}$ as a function of the targeted $c_{Y}(\tau)$:

[TABLE]

This permits to compute $h^{(m,\tau)}$ and $I^{(m,p,\tau)}$ analytically, using eqs.(12) and (13).

For stationary log-normal process $Y_{F_{2}}$, $F_{2}$ is not bijective. The entropy $H$ can be computed analytically, but this is not the case for $I^{(m,p,\tau)}$ and $h^{(m,\tau)}$.

III-B3 Correlation structures: short-range versus long range dependencies

To study the impact of the correlation structure, processes with short-range dependencies are compared against processes with long-range dependencies, the former being well-represented by an exponentially decaying covariance,

[TABLE]

and the latter by an algebraically decaying covariance,

[TABLE]

with $0.5<{\cal H}<1$.

III-C Two-point dependencies: Auto mutual information versus correlation

For both short-range and long range correlations, we compare the two LN processes obtained with $F_{1}$ and $F_{2}$, denoted $Y_{F_{1}}$, $Y_{F_{2}}$, against a Gaussian process with same correlation function, denoted $Y_{{\mathbb{I}}}$, where ${\mathbb{I}}$ stands for the identity transform. First order statistics (probability density function) and second order statistics measured by the correlation function or the power spectral density as well as two point statistics (i.e., $m=p=1$) entropy, entropy rate and AMI estimates are compared on Fig. 1 for $Y_{{\mathbb{I}}}$, $Y_{F_{1}}$ and $Y_{F_{2}}$ with short-range dependences and on Fig.2 for long-range dependences.

As expected, for each pair of short range and long range processes, the two log-normal processes, $Y_{F_{1}}$ and $Y_{F_{2}}$ cannot be distinguished one from the other using either PDF, correlation or power spectral estimation. Entropy being fully prescribed by the sole PDF, it is further not surprising that the two log-normal processes are also undistinguishable when using entropy estimates $H$ (Fig (1,d)).

However, and interestingly, the estimates $I^{(1,1,\tau)}$ and $h^{(1,\tau)}$ of two-point AMI and entropy rate show differences as functions of $\tau$ between processes $Y_{F_{1}}$ and $Y_{F_{2}}$, both for the short-range and long-range correlation cases. $I^{(1,1,\tau)}$ decreases and $h^{(1,\tau)}$ increases when $\tau\rightarrow\infty$, a predictable result as all dependencies vanish for large $\tau$ for stationary processes. Both functions (AMI and entropy rate) for processes $Y_{F_{1}}$ and $Y_{F_{2}}$ differ for small $\tau$, thus quantifying the difference between the joint distributions of these two processes, despite identical correlation functions and marginal distributions.

Fig. 1 and Fig. 2 show that the two-point AMI estimates $\hat{I}^{(1,1,\tau)}$ for log-normal process $Y_{F_{1}}$ ressembles far more to that of the Gaussian process $Y_{\mathbb{I}}$ than to that of the companion log-normal process $Y_{F_{2}}$, even they all three have the same correlation function. Very interestingly, it provides a deep insight into the temporal dynamics of the time series: processes $Y_{F_{1}}$ and $Y_{F_{2}}$ have same marginal distributions and same correlations functions, but the temporal dynamics of $Y_{F_{1}}$ is close to that of the Gaussian process $Y_{{\mathbb{I}}}$ while that of $Y_{F_{2}}$ is not. This also illustrates the benefits of AMI ${I}^{(1,1,\tau)}$ over entropy rate ${h}^{(1,\tau)}$ that mixes the impact of both marginal distribution and joint dynamical structure, thus making $Y_{F_{2}}$ comparable to $Y_{{\mathbb{I}}}$ at small $\tau$ when dependencies matter and $Y_{F_{2}}$ comparable to $Y_{F_{1}}$ at large $\tau$ when only the marginal distribution contribute to ${h}^{(1,\tau)}$.

Figs 1,d,e,f and 2,d,e,f also report theoretical $H$, $h$ and AMI for processes $Y_{\mathbb{I}}$ and $Y_{F_{1}}$, while theoretical $h$ and AMI are not available for process $Y_{F_{2}}$. All estimates well match the predicted theoretical behaviors (solid lines in Figs 1, 2). The theoretical behaviors of AMI and entropy rate of process $Y_{F_{1}}$ are obtained using eqs. (12) and (13) with the analytical expression (17).

These observations essentially hold both for the short-range (Fig. 1) and long-range (Fig. 2) correlations and yield the first key conclusion of this work: as opposed to correlation or spectral estimation, estimates of the two-points AMI (${I}^{(1,1,\tau)}$) and entropy rate (${h}^{(1,\tau)}$) are able to quantify differences in the temporal dynamics of processes even when they share the exact same marginal distributions and correlation functions, as long as they have different two-point joint statistics.

IV Estimation performance and illustrations

The aim of this section is to assess, by means of Monte Carlo simulations, the performance of the estimators $\hat{H}$, $\hat{I}^{(m,p,\tau)}$ and $\hat{h}^{(m,\tau)}$ for processes with different statistical structures and to study the impacts of key parameters either entering the definition of the estimators (number of nearest neighbors $k$) or characterizing the dependence structures of the processes (strength of the correlation $c_{Y}(\tau)$) as well as with sample sizes $N$.

IV-A Estimation performance

IV-A1 Monte Carlo simulations Set-up

The processes are synthetized numerically using circulant matrix embedding strategies, with the toolbox designed in [54, 53] and made available at www.hermir.org. As we have access to analytical curves only for $Y_{\mathbb{I}}$ and $Y_{F_{1}}$, these are the processes used to characterize the performance of the estimators.

For short and long memory, parameters are set $f_{\rm c}=0.1$ and ${\cal H}=0.7$ respectively. All variances are set to $\sigma_{Y}\equiv 1$. For each case, $100$ realizations of processes are synthesized and performances are assessed as averages (and standard deviations) across realizations.

To study the impact of the marginal distribution on estimation performance, process $Y_{F_{1}}$ is compared to process $Y_{\mathbb{I}}$.

IV-A2 Dependence with the sample size $N$

Fig. 3 reports performances as functions of the sample size $N$, for $2^{8}\leq N\leq 2^{16}$, for a fixed $k=5$. It clearly shows that the three estimators are asymptotically consistent (with vanishing biases and variances) irrespective of the marginal distribution and of the correlation being short-range or long-range. Unsurprisingly, variances are larger for the log-normal processes than they are for the Gaussian ones. Also, long-range dependencies imply much larger variances compared to short-range dependencies for a same $c_{Y}(\tau)=0.32$ (which corresponds to $\tau=1$ for the long range correlation process and $\tau=11$ for the short range one), thus clearly showing that long memory does not induce biases but does imply larger variances. The proposed estimators are overall robust to marginal distributions and dependence structures [46, 47], and have satisfactory performance as soon as $N\geq 2^{10}$, hence $N=2^{16}$ is used hereafter.

IV-A3 Dependence with the number of nearest neighbors $k$

Fig. 4 reports estimation performance as functions of $k$, for $N=2^{16}$ and same correlation $c_{Y}=0.32$ as above. Performances barely depend on $k$, but for a mild increase of the biases, and this irrespective of the marginal distributions and of short-range or long-range nature of the correlations. Because increasing $k$ further implies increasing the computational load, $k=5$ is used hereafter.

IV-A4 Influence of the correlation strength $c_{Y}(\tau)$

Fig. 5 reports estimation performance as a function of $c_{Y}(\tau)$ and shows that all biases slightly increase when $c_{Y}(\tau)$ increases. The variance appears to increase a little faster than linearly in $c_{Y}$. For both short-range and long-range correlations, the log-normal processes show larger standard deviations compared to the Gaussian ones, yet with similar biases. Further, for a given correlation level $c_{Y}(\tau)$ the long-range correlation structure shows larger standard deviations compared to the short-range correlation.

From this characterization, we can conclude that the entropy, AMI and entropy rate estimators show a low bias and standard deviation for the combination $k=5$ and $N=2^{16}$ of parameters used along the article. In addition, we can conclude that the standard deviation of all the three estimators is larger for strong correlations and non-Gaussian distributions, but their biases depend slightly on the marginal distribution and the correlation strength.

IV-B Process discrimination with auto mutual information

In Figs 1,e,f and 2,e,f, we show that both AMI and entropy rate across scales are able to discriminate between processes $Y_{F_{1}}$ and $Y_{F_{2}}$. The discrimination is far more effective for the smaller scales for both correlation structures (Long and short range).

To quantify the differences between $Y_{F_{1}}$ and $Y_{F_{2}}$ that we found in Figs 1 and 2, we first present the histograms of the estimations of the entropy rate (Fig. 6a) and of the AMI (Fig. 6b) for fixed $N=2^{16}$ and $\tau=1$. We perform a statistical ranksum test using Mann-Whitney U-test (Wilcoxon rank-sum test) and compute the $p$-value, for varying $\tau$ and $N$. As seen in Fig. 6d, the $p$-value quickly decreases with increasing $N$, as expected.

Fig. 6c quantifies the distinguishability between both log-normal processes across scales, showing that the discrimination is far better for smaller scales. This confirms the observations of Figs 1,e,f and 2,e,f.

IV-C Dependence structure beyond correlations

This sections further explores the ability of AMI ${I}^{(m,p,\tau)}$, with $m\geq 1$ and/or $p\geq 1$ to probe multiple-point dependences and hence fine details in temporal dynamics.

Fig. 7 reports higher-order AMI estimates $\hat{I}^{(m,p,\tau)}$ for several combination of $m$ and $p$, as functions of the embedding dimension $p$ for both short-range (left) and long-range (right) correlation functions, for processes $Y_{\mathbb{I}}$, $Y_{F_{1}}$ and $Y_{F_{2}}$. From top to bottom the $m$ value is equal to $1$, $2$, $3$ and $p$. Consequently, the first three lines study the behaviour of AMI on $p$ for fixed value of $m$, while the fourth line shows the behavior of AMI for $m=p$.

While two-point AMI yielded essentially identical conclusions for short-range and long-range correlation processes, higher order AMI reveals clear differences: For short-range correlation processes (left column), $Y_{\mathbb{I}}$, $Y_{F_{1}}$ and $Y_{F_{2}}$, $\hat{I}^{(m,p,\tau)}$ and ${I}^{(m,p,\tau)}$ when it can be computed analytically, show a very mild dependence of AMI with $m$ and $p$ for a fixed $\tau$. To the converse, for long-range correlation processes (right column), $Y_{\mathbb{I}}$, $Y_{F_{1}}$ and $Y_{F_{2}}$, $\hat{I}^{(m,p,\tau)}$ and ${I}^{(m,p,\tau)}$ when it can be computed analytically, show constant increase of AMI when $p$ or $m$ increase, for a fixed $\tau$. This fundamental observation clearly indicates that for short-range correlation, whatever $\tau$, there is little or no gain in the analysis of temporal dynamics brought by the study of higher order AMI. Conversely, for long-range correlation processes, ${I}^{(m,p,\tau)}$ captures more and more information in the temporal dynamics as the embedding dimensions are increased thus revealing a significant difference between short-range and long-range correlation structures. Analytical formulae for ${I}^{(m,p,\tau)}$ when $m,p>1$, are only available for the identity ${\mathbb{I}}$ and the natural transform $F_{1}$, namely for $Y_{{\mathbb{I}}}$ and $Y_{F_{1}}$.

In order to quantify the evolution of AMI when the embedding increases, we compute the $p$-value of Wilcoxon test between the statistic of $\hat{I}^{(1,1,1)}$ and $\hat{I}^{(m,p,1)}$ for the same processes as in Fig. 7. Fig. 8 quantifies the increase of AMI when the embedding increases. This evolution of AMI is far more important in the long range correlations case, as expected from Fig. 7. Further, table I shows that, for fixed $m+p$, AMI evolves differently depending on the precise values of $m$ and $p$ separately. For a fixed $m+p$ AMI is larger when $m$ and $p$ are close, which can be understood as a closer similarity between the two embedded vectors appearing as the arguments of the AMI in eq.(4).

V Conclusions and perspectives

We explored the evolution of AMI ($I^{(m,p,\tau)}$) and entropy rate ($h^{(m,\tau)}$) on the scale parameter $\tau$. Even when considering only two point interactions ($m=p=1$) AMI and entropy rate probe statistics of any order, and hence non-linear dependences, and as such appear as unambiguous generalizations of the correlation function and the power spectrum.

As an illustration, we analyzed two log-normal processes with identical 1-point distribution and identical correlation function. AMI ($I^{(1,1,\tau)}$) and entropy rate ($h^{(1,\tau)}$) are able to discriminate these two processes, which appear as identical for classical analysis, by enlightening differences in high order dependences. The behavior of these quantities along scales reveal the existence of stronger high-order dependences at smaller scale, which allows an easier separation of processes.

In addition, AMI generalizes easily to consider explicitly interactions between more than two points. We thus probe the complexity of the dependence structure above linear order, i.e., the additional information measured by AMI when considering an extra point should indicate the existence and relevance of next-order interactions. As an example, comparing the effect of increasing the order of interactions for two dependence structures showed that while exponential decay — short range — dependences do not involve next order interactions, power law — long range — dependences do.

The same qualitative results were obtained when considering others one-point distribution, which supports the validity of the methodology for any distribution and dependence structure. This generality should make the methodology very interesting to perform non-Gaussian processes characterization in several different applications.