Spike Spectra for Recurrences

Abstract

In recurrence analysis, the τ-recurrence rate encodes the periods of the cycles of the underlying high-dimensional time series. It, thus, plays a similar role to the autocorrelation for scalar time-series in encoding temporal correlations. However, its Fourier decomposition does not have a clean interpretation. Thus, there is no satisfactory analogue to the power spectrum in recurrence analysis. We introduce a novel method to decompose the τ-recurrence rate using an over-complete basis of Dirac combs together with sparsity regularization. We show that this decomposition, the inter-spike spectrum, naturally provides an analogue to the power spectrum for recurrence analysis in the sense that it reveals the dominant periodicities of the underlying time series. We show that the inter-spike spectrum correctly identifies patterns and transitions in the underlying system in a wide variety of examples and is robust to measurement noise.

Keywords: bifurcations; decomposition; frequency analysis; recurrence analysis.

Figure A1

Subset of the frequency time series of (A) Central Europe (CE) [45] and (C) Great Britain (GB) [44] along with their autocorrelation functions in (B,D), respectively. The normalized autocorrelations in (B,D) were computed with the entire time series, whereas the subsets shown in (A,C) contain only 1441 samples.

Figure A2

Fourier spectra of the frequency entire time series from (A) Central Europe (CE), $N_{\mathrm{CE}}=587,520$ and (B) Great Britain (GB), $N_{\mathrm{GB}}=1,581,120$. A subset of the underlying time series is shown in Figure A1.

Figure A3

Same as in Figure 4, but here with 5% additive Gaussian white noise on each component x, y and z. (A) Trajectory of the system in a period-2 (parameter $a=0.36$), (B) in a period-3 (parameter $a=0.41$) and (C) in a chaotic regime (parameter $a=0.428$). (D–F) The corresponding RPs, obtained by using a recurrence threshold corresponding to a 10% global recurrence rate for (D,E) and 5% for (F). (G–I) $\tau$-RRs of the shown RPs, (J–L) the corresponding inter-spike spectra, and (M–O) the Fourier power spectra. The appearance of the inter-spike spectra in (J–L), and the Fourier spectra in (M–O) are unaffected by the additive noise.

Figure A3

Same as in Figure 4, but here with 5% additive Gaussian white noise on each component x, y and z. (A) Trajectory of the system in a period-2 (parameter $a=0.36$), (B) in a period-3 (parameter $a=0.41$) and (C) in a chaotic regime (parameter $a=0.428$). (D–F) The corresponding RPs, obtained by using a recurrence threshold corresponding to a 10% global recurrence rate for (D,E) and 5% for (F). (G–I) $\tau$-RRs of the shown RPs, (J–L) the corresponding inter-spike spectra, and (M–O) the Fourier power spectra. The appearance of the inter-spike spectra in (J–L), and the Fourier spectra in (M–O) are unaffected by the additive noise.

Figure A4

Wasserstein distances for the obtained inter-spike spectra of the $\tau$-RRs for additive noise levels of up to 50% for the discussed Rössler dynamics, (A) period-2 limit-cycle, (B) period-3 limit-cycle and (C) chaos. For each noise level the obtained spectrum is compared to the noise-free spectrum (these are shown in Figure 4J–L for a regularization threshold corresponding to $\rho =0.95$) by computing its Wasserstein distance. Inter-spike spectra were obtained with a LASSO regression and three different regularization thresholds corresponding to $\rho =0.9$ (blue), $\rho =0.95$ (orange) and $\rho =0.99$ (green) accordance of $\tau$-RRs and re-composed signals. Up to a noise level of 20% the distances are small and rather constant. The fluctuations depend on the chosen regularization threshold as well as on the underlying dynamics.

Figure A5

Pearson correlation coefficient ${\rho }_{\overrightarrow{s},\tilde{\overrightarrow{s}}}$ between the input signal signal $\overrightarrow{s}$ and the re-composed signal $\tilde{\overrightarrow{s}}={\mathbf{X}}^T\widehat{\mathbf{\beta }}$ as a function of the regularization parameter $\alpha$ for (A) LASSO and (B) sequentially thresholded least squares (STLS) regression (see also Section 2). The input time series is of length $N=200$ and stems from the $\tau$-RR of Rössler system in regular dynamics.

Figure 1

Schematic illustration of a $\tau$-recurrence-rate-based spectrum. (A) x-component time series of the Lorenz63-System (Equation (A1)) and (B) its corresponding Fourier power spectrum. (C) Reconstructed state space portrait from the time series shown in (A) using PECUZAL time-delay embedding [16]. (D) Subset of the recurrence plot and the corresponding $\tau$-recurrence rate obtained from the state space trajectory in (C). (E) Fourier Power spectrum obtained from the $\tau$-recurrence rate (subset shown in panel (D)) [15]. (D,E) show the results of a part of the time series, which is highlighted in pink in (A).

Figure 2

The transformation of a Dirac comb (series of Dirac delta functions) with a single inter-spike period $T_{\mathrm{is}}=100$ ($\widehat{=}f=0.01$) into the frequency domain. (A) Dirac Comb (DC) with equal amplitudes and (B) its FFT-based power spectral density. (C) Proposed inter-spike spectrum of the signal in (A) showing a single frequency, which corresponds to the inter-spike period $T_{\mathrm{is}}$ ($f=0.01$). (D) DC with randomly chosen amplitudes and same $T_{\mathrm{is}}$ as in (A), and (E) is its FFT-based power spectral density. (F) Proposed inter-spike spectrum of the signal in (D) showing a single frequency, which corresponds to the expected inter-spike period $T_{\mathrm{is}}$ ($f=0.01$). Inter-spike spectra were obtained with a LASSO regression and a regularization threshold corresponding to $\rho =0.9$ accordance of the signals in (A,D) and its re-composed signals (c.f. Section 2).

Figure 3

(A) Example of a full set of basis functions for an input signal of length $N=4$, aligned in the matrix ${\mathbf{X}}^{\prime }$. Inter-spike periods $T_{\mathrm{is}}$ larger than $\lceil N/2\rceil +1$ lead to redundant basis functions (i.e., repeated lines in X, red sheared) or basis functions, which cannot be uniquely assigned to a certain inter-spike period (blue sheared). $T_{\mathrm{is}}$ for each row can be obtained by Equation (6). (B) The final set of unambiguous, but still linearly dependent, basis functions aligned in the matrix X.

Figure 4

Inter-spike spectra of the $\tau$-RR of the Rössler system in three different dynamical regimes with parameters $b=2$, $c=4$. (A) Trajectory of the system in a period-2 (parameter $a=0.36$), (B) in a period-3 (parameter $a=0.41$) and (C) in a chaotic regime (parameter $a=0.428$). (D–F) The corresponding RPs, obtained by using a recurrence threshold corresponding to a 10% global recurrence rate for (D,E) and 5% for (F). (G–I) $\tau$-RRs of the shown RPs. (J–L) The proposed inter-spike spectra of the $\tau$-RRs shown in panels (G–I). Spectra were obtained with a LASSO regression and a regularization threshold corresponding to $\rho =0.95$ accordance of $\tau$-RRs and re-composed signals. The distance ratio of the peaks reflect the limit cycle dynamic. (M–O) Fourier power spectra of the $\tau$-RRs shown in panels (G–I).

Figure 5

Peak positions of the obtained inter-spike spectra of the $\tau$-RRs for additive noise levels up to 50% for the discussed Rössler dynamics, (A) period-2 limit-cycle, (B) period-3 limit-cycle and (C) chaos. The size of the plotted markers scale with the detected peak height. The noise-free spectra are shown in Figure 4J–L and an example of these spectra with 5% additive noise is shown in Figure A3J–L. Spectra were obtained with a LASSO regression and a regularization threshold corresponding to $\rho =0.95$ accordance of $\tau$-RRs and re-composed signals.

Figure 6

(A) Bifurcation diagram and Lyapunov exponent and of the Logistic map as a function of the control parameter r. (B) Number of significant peaks ($\alpha =0.05$) in the inter-spike spectrum of the $\tau$-RR and its Pearson correlation coefficient to the Lyapunov exponent shown in (A) (white noise surrogates). (C) Same as (B), but for iterative Amplitude Adjusted Fourier Transform (iAAFT) surrogates [26,27]. To obtain the inter-spike spectra, we used a LASSO regression and a regularization threshold corresponding to the $\rho =0.95$ accordance of $\tau$-RRs and re-composed signals.

Figure 7

Averaged evolutionary Fourier power spectra of recorded power grid frequency time series of (A) Central Europe (CE) and (B) Great Britain (GB) (see Figure A1 in Appendix B). The corresponding averaged inter-spike spectra of the according $\tau$-RRs, Equation (2), are shown in (C,D), respectively. Vertical red dashed lines correspond to 15, 30 and 60 $\min$. For technical details on the calculation of the spectra shown, i.e, the preprocessing and window size being used, the reader is referred to Appendix B.

Figure 8

(A) Evolutionary Fourier power spectra of eccentricity time series. (B) Inter-spike spectrogram of the $\tau$-recurrence rate of the eccentricity time series and (C) its Fourier spectrogram. Horizontal black dashed lines highlight the analytically expected orbital periods of 405, 124 and 95 kyrs. For comparability, in all cases the spectra aligned in the columns of the shown plots are normalized to probabilities (sum of unity for each power spectrum). For further computational details, please refer to the main text.