acf
- mdtools.statistics.acf(x, axis=None, dt=1, dtau=1, tau_max=None, center=True, unbiased=False)[source]
Calculate the autocorrelation function of an array.
- Parameters:
x (
array_like
) – The input array.axis (
int
orNone
, optional) – The axis along which to calculate the ACF. By default, the flattened input array is used.dt (
int
orNone
, optional) – Distance between restarting times (actually restarting indices) \(t_0\). IfNone
is given, dt is set to the current lag time \(\tau\) which results in independent sampling windows.dtau (
int
, optional) – Distance between evaluated lag times (actually lag indices) \(\tau\).tau_max (
int
orNone
, optional) – Maximum lag time \(\tau\) upon which to compute the ACF. Because the first lag time is 0, the last evaluated lag time is actuallytau_max - 1
. By default, tau_max is set tox.size
is axis isNone
and tox.shape[axis]
otherwise. tau_max must lie within the interval[0, len(x)]
.center (
bool
, optional) – IfTrue
, center the input array around its mean, i.e. subtract the sample mean \(\bar{X}\) when calculating the variance and covariance.unbiased (
bool
, optional) – IfTrue
, the covariance \(\text{Cov}[X_{t_0}, X_{t_0 + \tau}]\) is normed by the actual number of sample points \(t_0\) (which depends on the lag time \(\tau\)). IfFalse
, the covariance is for all lag times normed by the number of sampling points at \(t_0 = 0\) (which is equal to the length of the input array x). Note that setting unbiased toTrue
might result in values of the ACF that lie outside the interval [-1, 1], especially for lag times larger thanlen(x) // 2
if center isTrue
.
- Returns:
ac (
numpy.ndarray
) – Autocorrelation function of x. If axis isNone
, the output array has shape(n,)
wheren
isint(np.ceil(tau_max / dtau))
. If a specific axis is given, ac has the same shape as x except along the given axis where the length is given byn
.- Raises:
ValueError : – If the first element of ac is not unity or one or more elements of ac fall outside the interval [-1, 1].
See also
mdtools.statistics.acf_np()
Different implementation of the ACF using
numpy.correlate()
mdtools.statistics.acf_se()
Calculate the standard errors of an autocorrelation function
mdtools.statistics.acf_confint()
Calculate the confidence intervals of an autocorrelation function
Notes
This function computes the autocorrelation function (ACF) of the input array x assuming that x describes a weakly stationary stochastic process. That means, the mean does not vary with time, the variance is finite for all times and the autocovariance does not depend on the starting time \(t_0\) but only on the lag time \(\tau\). In this case the ACF is given by:[1],[2]
\[C(\tau) = \frac{\text{Cov}[X_{t_0}, X_{t_0 + \tau}]}{\text{Var}[X]} = \frac{\langle (X_{t_0} - \mu) \cdot (X_{t_0 + \tau} - \mu) \rangle}{\langle (X - \mu)^2 \rangle}\]The ACF strictly fulfills \(C(\tau) \in [-1, 1] \text{ } \forall \text{ } \tau \in \mathbb{R}\). Numerically, the ACF is calculated using the following formula:[3],[4]
\[C_\tau = \frac{\frac{1}{T^*} \sum_{t_0=0}^{T-\tau} (X_{t_0} - \bar{X}) (X_{t_0 + \tau} - \bar{X})}{\frac{1}{T} \sum_{t_0=0}^T (X_{t_0} - \bar{X})^2} = \frac{\frac{1}{T^*} \sum_{t_0=0}^{T-\tau} (X_{t_0} - \bar{X}) (X_{t_0 + \tau} - \bar{X})}{C_0}\]The increment of the sums in the formula of \(C_\tau\) is given by dt. \(\bar{X} = \frac{1}{T} \sum_{t=1}^T X_t\) is the sample mean. For the sample mean, the increment of the sum is always 1, independent of the choice of dt. Note that the sample mean is only subtracted when center is
True
. T is always the length of the input array x, independent of the choice of tau_max. \(T^*\) is either \(T - \tau\) if unbiased isTrue
or simply \(T\) if unbiased isFalse
.For the sake of completeness, the ACF of a non-stationary stochastic process is given by:[1],[2]
\[C(t_1, t_2) = \frac{\text{Cov}[X_{t_1}, X_{t_2}]}{\sqrt{\text{Var}[X_{t_1}]} \cdot \sqrt{\text{Var}[X_{t_2}]}} = \frac{\langle (X_{t_1} - \mu_{t_1}) \cdot (X_{t_2} - \mu_{t_2}) \rangle}{\sqrt{\langle (X_{t_1} - \mu_{t_1})^2 \rangle} \cdot \sqrt{\langle (X_{t_2} - \mu_{t_2})^2 \rangle}}\]References
Examples
>>> a = np.array([-2, 2, -2, 2, -2]) >>> mdt.stats.acf(a) array([ 1. , -0.8 , 0.56666667, -0.4 , 0.13333333]) >>> mdt.stats.acf(a, dt=None) array([ 1. , -0.8 , 0.26666667, -0.2 , 0.13333333]) >>> mdt.stats.acf(a, unbiased=True) array([ 1. , -1. , 0.94444444, -1. , 0.66666667]) >>> mdt.stats.acf(a, unbiased=True, dt=None) array([ 1. , -1. , 0.66666667, -1. , 0.66666667]) >>> mdt.stats.acf(a, center=False) array([ 1. , -0.8, 0.6, -0.4, 0.2]) >>> mdt.stats.acf(a, center=False, dtau=2) array([1. , 0.6, 0.2]) >>> mdt.stats.acf(a, center=False, tau_max=len(a)//2) array([ 1. , -0.8]) >>> mdt.stats.acf(a, center=False, dt=None) array([ 1. , -0.8, 0.4, -0.2, 0.2]) >>> mdt.stats.acf(a, center=False, unbiased=True) array([ 1., -1., 1., -1., 1.]) >>> mdt.stats.acf(a, center=False, unbiased=True, dt=None) array([ 1., -1., 1., -1., 1.])
>>> a = np.array([-2, -2, -2, 2, 2]) >>> mdt.stats.acf(a) array([ 1. , 0.36666667, -0.26666667, -0.4 , -0.2 ]) >>> mdt.stats.acf(a, dt=None) array([ 1. , 0.36666667, -0.06666667, -0.2 , -0.2 ]) >>> mdt.stats.acf(a, unbiased=True) array([ 1. , 0.45833333, -0.44444444, -1. , -1. ]) >>> mdt.stats.acf(a, unbiased=True, dt=None) array([ 1. , 0.45833333, -0.16666667, -1. , -1. ]) >>> mdt.stats.acf(a, center=False) array([ 1. , 0.4, -0.2, -0.4, -0.2]) >>> mdt.stats.acf(a, center=False, dt=None) array([ 1. , 0.4, 0. , -0.2, -0.2]) >>> mdt.stats.acf(a, center=False, unbiased=True) array([ 1. , 0.5 , -0.33333333, -1. , -1. ]) >>> mdt.stats.acf(a, center=False, unbiased=True, dt=None) array([ 1. , 0.5, 0. , -1. , -1. ])
>>> a = np.array([[-2, 2, -2, 2, -2], ... [-2, 2, -2, 2, -2], ... [-2, -2, -2, 2, 2], ... [-2, -2, -2, 2, 2]]) >>> mdt.stats.acf(a, center=False) array([ 1. , -0.15, 0. , -0.25, -0.2 , 0.55, 0. , 0.15, -0.1 , -0.25, 0.1 , 0.05, 0.1 , 0.05, -0.2 , 0.05, 0. , 0.05, 0. , -0.05]) >>> mdt.stats.acf(a, center=False, axis=0) array([[ 1. , 1. , 1. , 1. , 1. ], [ 0.75, 0.25, 0.75, 0.75, 0.25], [ 0.5 , -0.5 , 0.5 , 0.5 , -0.5 ], [ 0.25, -0.25, 0.25, 0.25, -0.25]]) >>> mdt.stats.acf(a, center=False, axis=0, unbiased=True) array([[ 1. , 1. , 1. , 1. , 1. ], [ 1. , 0.33333333, 1. , 1. , 0.33333333], [ 1. , -1. , 1. , 1. , -1. ], [ 1. , -1. , 1. , 1. , -1. ]]) >>> mdt.stats.acf(a, center=False, axis=1) array([[ 1. , -0.8, 0.6, -0.4, 0.2], [ 1. , -0.8, 0.6, -0.4, 0.2], [ 1. , 0.4, -0.2, -0.4, -0.2], [ 1. , 0.4, -0.2, -0.4, -0.2]]) >>> mdt.stats.acf(a, center=False, axis=1, unbiased=True) array([[ 1. , -1. , 1. , -1. , 1. ], [ 1. , -1. , 1. , -1. , 1. ], [ 1. , 0.5 , -0.33333333, -1. , -1. ], [ 1. , 0.5 , -0.33333333, -1. , -1. ]])
>>> a = np.array([[[-2, 2, -2, 2, -2], ... [-2, 2, -2, 2, -2]], ... ... [[-2, -2, -2, 2, 2], ... [-2, -2, -2, 2, 2]]]) >>> mdt.stats.acf(a, center=False) array([ 1. , -0.15, 0. , -0.25, -0.2 , 0.55, 0. , 0.15, -0.1 , -0.25, 0.1 , 0.05, 0.1 , 0.05, -0.2 , 0.05, 0. , 0.05, 0. , -0.05]) >>> mdt.stats.acf(a, center=False, axis=0) array([[[ 1. , 1. , 1. , 1. , 1. ], [ 1. , 1. , 1. , 1. , 1. ]], [[ 0.5, -0.5, 0.5, 0.5, -0.5], [ 0.5, -0.5, 0.5, 0.5, -0.5]]]) >>> mdt.stats.acf(a, center=False, axis=0, unbiased=True) array([[[ 1., 1., 1., 1., 1.], [ 1., 1., 1., 1., 1.]], [[ 1., -1., 1., 1., -1.], [ 1., -1., 1., 1., -1.]]])
>>> mdt.stats.acf(a, center=False, axis=1) array([[[1. , 1. , 1. , 1. , 1. ], [0.5, 0.5, 0.5, 0.5, 0.5]], [[1. , 1. , 1. , 1. , 1. ], [0.5, 0.5, 0.5, 0.5, 0.5]]]) >>> mdt.stats.acf(a, center=False, axis=1, unbiased=True) array([[[1., 1., 1., 1., 1.], [1., 1., 1., 1., 1.]], [[1., 1., 1., 1., 1.], [1., 1., 1., 1., 1.]]]) >>> mdt.stats.acf(a, center=False, axis=2) array([[[ 1. , -0.8, 0.6, -0.4, 0.2], [ 1. , -0.8, 0.6, -0.4, 0.2]], [[ 1. , 0.4, -0.2, -0.4, -0.2], [ 1. , 0.4, -0.2, -0.4, -0.2]]]) >>> mdt.stats.acf(a, center=False, axis=2, unbiased=True) array([[[ 1. , -1. , 1. , -1. , 1. ], [ 1. , -1. , 1. , -1. , 1. ]], [[ 1. , 0.5 , -0.33333333, -1. , -1. ], [ 1. , 0.5 , -0.33333333, -1. , -1. ]]])