There are nonequilibrium thermodynamic systems that for short time periods are in equilibrium states, but their parameters vary over longer time scales. An elementary example of such a system is a Brownian particle moving through a slowly fluctuating environment. If the Brownian particle is always in local equilibrium, then its velocity has a Gaussian distribution with slowly varying parameters. The average global velocity distribution is a superposition of Gaussians weighted by the probability density of the slowly fluctuating parameters. Depending on the type of the statistical distribution of the slow fluctuations of the environment, the particle velocity may have different types of statistical distributions named superstatistics beck03 .

A time series obtained by measuring a thermodynamic system with superstatistics can be modeled as the product of two stochastic processes

where $t=0,\pm 1,\pm 2,\mathrm{\dots }$ and $V_t$ and $Z_t$ are independent. The Gaussian white noise $\{Z_t\}$ models the equilibrium thermodynamic fluctuations, while the stochastic process $\{V_t\}$ with strictly positive values describes the slow fluctuations of the environment. If $\{Z_t\}$ is not only uncorrelated, but independent and identically distributed (i.i.d.), then $\{X_t\}$ is called heteroskedastic. We emphasize that the independence of $V_t$ and $Z_t$ does not preclude the dependence of $V_t$ on the preceding terms $Z_{t-s}$ and $V_{t-s}$, $s>0$, as for example, in the ARCH and GARCH models engle82 ; boller86 .

Although the financial markets are not usual thermodynamic systems, many of the methods in statistical physics can be applied to such complex social phenomena voit05 . For a financial time series the white noise $\{Z_t\}$ in the heteroskedastic process (1) models the efficiency of the financial markets, i.e., the uncorrelation of successive price variations. The significant long-range correlations of the absolute values of the price fluctuations are due to the slowly varying stochastic process $\{V_t\}$, which modulates the amplitude of the white noise fluctuations. As in finance, we name it volatility.

A numerical method to estimate the two factors in Eq. (1) is the partitioning of the time series into slices for which the thermodynamic local equilibrium holds straet09 . The partitioning is achieved by imposing the condition that the mean kurtoses of $\{X_t\}$ over the slices equal the Gaussian kurtosis. The superstatistics method has been applied to a variety of complex systems with time scale separation in turbulence, share price fluctuations, cosmic rays, traffic delays, metastasis and cancer survival, etc. beck09 . Because the normality is tested by means of the fourth moment, the superstatistics method is sensitive to the presence of outliers. For example, it was necessary to eliminate the extreme events from some financial time series in order to analyze them straet09 . A similar approach for intraday returns is presented in gerig09 . The main difference is that the slices length is predetermined to a day and the volatility is considered constant over intraday time scales.

The equilibrium thermodynamic fluctuations are not only Gaussian, but also independent, and, consequently, uncorrelated. Hence the separation of the components of the stochastic process (1) can be obtained by imposing instead of the normality condition, the condition that the absolute values of the estimated noise be uncorrelated. The question is about the relation between these two conditions, i.e., if we impose one condition, then to what degree is the other condition also satisfied. In this Rapid Communication we show that by imposing the uncorrelation condition, the estimated white noise becomes almost Gaussian, proving that it entails the normality condition.

Our results indicate that the uncorrelation condition is numerically more efficient than the normality one. We apply this method to the time series of the daily S&P500 index analyzed in straet09 , without its partitioning and without eliminating any extreme value. It is possible that this greater efficiency be due to a more complete exploitation of the temporal variation of all the observed values, not only the variation from one slice to the other.

Although the physical systems modeled by Eq. (1) are normal in the conventional thermodynamic sense, they seem anomalous by featuring non-Gaussian statistics. If the volatility $\{V_t\}$ is stationary, then $\{X_t\}$ is also stationary and its probability density function (pdf) is a scale mixture of normal distributions given by the integral

where $f_V$ is the pdf of $V_t$ andrews74 . According to this formula, $f_X(x)$ is the result of the superposition of Gaussians whose standard deviations $v$ are realizations of the volatility $V_t$, that is why we say that $f_X(x)$ is the mixture of Gaussians with scales distributed according to $f_V$. A large class of systems with polydisperse features has the pdf with power-law tails as a result of the superposition of Gaussians with the variances obeying a gamma distribution gherghiu04 .

A direct generalization of the scale mixture of normal distributions is obtained if the volatility is nonstationary. Since $f_V$ in Eq. (2) explicitly depends on time, then $f_X$ also varies in time. The global average pdf is the time average of $f_X(x,t)$

and we refer to it as a scale mixture of Gaussian white noises. This distinct denomination emphasizes the importance of the time dependence of the volatility which allows for a greater range of shapes for the pdfs than Eq. (2). In this case the volatility can describe a deterministic variation imposed by the environment on the measured system or a superposition of a stationary random fluctuations on such a deterministic evolution.

We determine the white noise $\{Z_t\}$ by first estimating the volatility $\{V_t\}$ by one of the two main methods in quantitative finance poon05 . The first method expresses the volatility by the explicit relations that define stochastic models of the financial time series. The second method is based on the averaging of historical values of the absolute or square returns at a high enough sampling frequency. Such estimators are called “realized volatility” and it is shown that for continuous stochastic processes they are unbiased estimators andersen02 . In practice various types of averaging are used. The best known is the exponential weighted moving average of the Riskmetrics model riskmetr96 . The daily volatility is also often estimated by averaging the intraday data andersen98 . In this paper we estimate the volatility by means of a modified form of the moving average presented in liu99 and we find the optimum moving average using the uncorrelation condition.

We consider a time series $\{x_t,t=1,2,\mathrm{\dots },T\}$ as a finite realization of the stochastic process $\{X_t\}$. We follow the practice of denoting by small letters the realizations of the stochastic processes denoted by the corresponding capitals. In this way we keep clear the distinction between the different nature of a stochastic process (family of measurable functions defined on a probability space) and one of its realizations which is a series of real numbers equal with the values taken by the measurable functions at particular events. In fact the stochastic process does not model the time series, but the random phenomenon which is assumed to generate the observed time series.

If $K$ is the semilength of the averaging window, then for we define the moving average

If $t⩽K$ ($t>T-K$), then the average is taken over the first $t+K$ (the last $T-t+K+1$) values of $|x_t|$. This asymmetric average forces the values near the time series boundaries to follow the variations of the interior values. If we consider ${\vartheta }_t(K)$ as a volatility estimator, then from Eq. (1) it follows that the estimator of the white noise is

We have to find the optimum value $K_0$ satisfying the condition that the series $\{|{\zeta }_t(K)|\}$ be uncorrelated. This condition cannot be imposed to the signed series $\{{\zeta }_t(K)\}$ since, like the initial time series $\{x_t\}$, it is uncorrelated. This property results from the randomness of the sign variations of the white noise $\{Z_t\}$ in Eq. (1), which is preserved by multiplication with the strictly positive values of the slowly varying volatility $\{V_t\}$.

The sample autocorrelation function (ACF) of an uncorrelated infinite time series is identical zero. Bartlett’s formula states that for a finite i.i.d. Gaussian process, the pdf of the sample ACF is Gaussian with mean $0$ and variance $1/T$, where $T$ is the length of the time series brockdav96 . The proof presented in brockdav96 can be reproduced with minimal changes for the absolute values of an i.i.d. Gaussian process and it follows that Bartlett’s formula holds in this case too.

We measure the deviation from normality of the sample ACF of the time series $\{\left|{\zeta }_t(K)\right|\}$, denoted ${\widehat{\rho }}_{\left|\zeta \right|}$, by the statistic used in the Kolmogorov-Smirnov test

where $\widehat{F}$ is the sample cumulative distribution function (cdf) of ${\widehat{\rho }}_{\left|\zeta \right|}$ and $G$ is the theoretical cdf of the normal distribution with Bartlett’s parameters. The quantity $\epsilon$ is an index measuring the nonnormality of ${\widehat{\rho }}_{\left|\zeta \right|}$, i.e., the serial correlation of $\{\left|{\zeta }_t(K)\right|\}$. If $K_0$ is the value for which $\epsilon$ is minimum, then the volatility estimator is ${\widehat{v}}_t={\vartheta }_t(K_0)$ and from Eq. (4) it follows that the estimator of the white noise is ${\widehat{z}}_t=x_t/{\widehat{v}}_t$.

We test this algorithm on the daily S&P500 index studied in straet09 on which several scaling analyses have been performed liu99 ; Gopi99 ; mantegna95 . It has also been proved by a statistical analysis that this time series can be modeled as the product of the volatility and an i.i.d. stochastic process starica05 . This conclusion is in accordance with the heteroskedastic decomposition in Eq. (1) on which all the volatility stochastic models are build poon05 . Its validity is sustained also by the results analogous to those for the S&P500 index obtained with the uncorrelation condition for several stock market indices (Dow Jones, Nikkei, FTSE, etc.) and several large cap stocks (INTC, MSFT, IBM, MRK, etc.).

We examine the S&P500 index in the interval between January 1st, 1950 and December 31st, 2009 containing $T=15097$ trading days. If we denote by $\{p_t,t=1,2,\mathrm{\dots },T\}$ the daily closing values of the S&P500 index, then the logreturn over the temporal interval $[t,t+\mathrm{\Delta }t]$ is defined as

For $\mathrm{\Delta }t>1$ the logreturns are computed only for $t=n\mathrm{\Delta }t$, otherwise the temporal intervals would overlap each other and spurious correlations would occur between the logreturns. Also, the mean is extracted from the final time series $\{r_t(\mathrm{\Delta }t)\}$.

We suppose that, for a given $\mathrm{\Delta }t$, the time series $\{r_t(\mathrm{\Delta }t),t=1,2,\mathrm{\dots },T\}$ is a finite sample of an infinite discrete stochastic process $\{R_t(\mathrm{\Delta }t),t=0,\pm 1,\pm 2,\mathrm{\dots }\}$ of the form of Eq. (1). The average (3) is performed only on the daily logreturns ($\mathrm{\Delta }t=1$), so that the information at the highest available frequency is used even if we estimate the volatility for $\mathrm{\Delta }t>1$. This choice is in agreement with the theoretical result that the realized volatility is closer to the real one if data frequency is increased andersen02 .

The index $\epsilon$ defined by Eq. (5) depends on two parameters ($K$ and $H$). The choice of the number $H$ of sample ACF values is a difficult problem with no simple solution escan09 . As a general rule $H$ should be smaller than a quarter of the time series length box94 . However, the sample ACF ${\widehat{\rho }}_{\left|\zeta \right|}$ is different from that of a white noise especially for the small values of $h$. Moreover, if the serial correlation for small $h$ is reduced, then the entire ACF becomes negligible. As a consequence, we limit the value of $H$ in Eq. (5) only to several tens.

In Fig. 1(a) we present the variation of the index $\epsilon$ for $\mathrm{\Delta }t=1$ with respect to the semilength of the averaging window for several values of $H$. Each curve has two minima. The first minimum corresponds to an averaging with $K=1$, which does not significantly damp the fluctuations of the logreturns, so that the average ${\vartheta }_t(K)$ preserves these fluctuations. We are interested in the second minimum occurring when $K>1$ for which the average ${\vartheta }_t(K)$ is slowly varying. For other values of $\mathrm{\Delta }t$ the first minimum is greater than the second one or even does not exist and then the minimum supplying the estimated volatility coincides with the global minimum.

The minimum values of the index $\epsilon$ for which the white noise is estimated are compared in Fig. 1(b) with ${\epsilon }_0$ which is the index $\epsilon$ computed for the initial logreturns time series. The minimum of $\epsilon$ is obtained by exhaustive search for $K⩽150$ and $H\in \{10,20,30,40,50\}$. One notices that for all time scales $\mathrm{\Delta }t$, ${\epsilon }_0$ is close to 1, indicating that all the $H$ values of the sample ACF implied in the computation of ${\epsilon }_0$ lay outside the variation range of the ACF of a white noise. The values of ${\epsilon }_{\min }$ are significantly smaller than the corresponding ${\epsilon }_0$ showing that the estimated white noise is much closer to an uncorrelated time series than the initial logreturns.

In order to verify that the uncorrelation condition entails the normality of the estimated white noise, we compute the Kolmogorov-Smirnov statistic $D$ of $\{{\widehat{z}}_t\}$. In comparison with Eq. (5) the difference is that $\widehat{F}$ is the cdf of the normalized estimated white noise, not of its sample ACF, and $G$ is the cdf of a normalized Gaussian. Figure 2 shows the values of $D$ for $\mathrm{\Delta }t⩽10$ compared with those for the normalized initial logreturns $D_0$. It also contains the mean and the standard deviation of $D$ for statistical ensembles of 1000 numerically generated i.i.d. Gaussian time series having the same length as $r_t(\mathrm{\Delta }t)$ of S&P500. The results show that the estimated white noise $\{{\widehat{z}}_t(\mathrm{\Delta }t)\}$ has a probability distribution much closer to a Gaussian than that of the initial logreturns. One notices that for $\mathrm{\Delta }t=2$ it cannot be differentiated from an i.i.d. Gaussian time series.

Figure 3 shows the logreturns, the estimated volatility, and the estimated white noise after the year 2000 for $\mathrm{\Delta }t=2$. The increased variability of the logreturns amplitude (volatility clustering) is no more present in the estimated white noise. Even the most volatile period of the financial crisis which started in 2008 is correctly described by the estimated volatility. We have obtained this result without removing any extreme value from the logreturns as it is necessary in the superstatistics method straet09 .

In order to test the separation algorithm in its simplest form, we have processed the entire logreturns time series, without any partitioning. However, for such a long time series ($T=15097$) it is possible that the search for the minimum of the index $\epsilon$ be more efficient over shorter fragments. Then the white noise obtained by joining the white noises estimated over shorter segments could be more similar to an uncorrelated time series. The partitioning of a nonstationary time series might be obtained by existing algorithms Fukuda04 .

The method presented in this paper numerically separates the two components of a scale mixture of Gaussian white noises and it can be applied to many phenomena characterized by the modulation in amplitude of the thermodynamic equilibrium fluctuations by a slowly varying process. Such phenomena are those already studied by means of the superstatistics method: turbulence, share price fluctuations, cosmic rays, traffic delays, metastasis and cancer survival, etc. beck09 . Some biophysical processes, for example, the human heart rate fluctuations, also present these features Galvan01 .

If the white noise $\{Z_t\}$ does not have a Gaussian distribution, Bartlett’s formula has to be modified for each type of pdf. For instance, at time scales of minutes the logreturns have a two-tailed exponential distribution Silva04 . In such cases the separation of the components using the uncorrelation condition gives different results depending on the assumed type of the probability distribution of the white noise. The choice of the correct pdf might be possible by means of the resemblance degree between the pdf of the estimated noise and the one that was initially assumed.