Recurrence plots

Written by Paul Bourke
January 1998


Recurrence plots are used to reveal non-stationarity of a series as well as to indicate the degree of aperiodicity. They were first proposed by Eckmann et al around 1986 in order to study phase space orbits. For example, for a stationary system the recurrence plot should generally be homogeneous along the diagonal.

Consider a series x with N terms

x0, x1, x2, x3, x4, x5, .... xi, ..... xN-1

Form all the vectors yi of dimension (length) D with lag (delay) d.

yi = ( xi, xi+d, xi+2d, .... )

D >= 2 and d >= 1

This is normally refered to as embedding x in dimension D with lag d.

The recurrence plot is formed by comparing all embedded vectors with each other, and drawing points when the distance between two vectors is below some threshold. More precisely, we draw a point at coordinate (i,j) if the i'th and j'th embedded vectors are less than some distance r apart, eg: a point is drawn if

||yi - yj|| < r

i is plotted along the horizontal axis, j on the vertical axis.

If the recurrence plot contains a homogenous but irregular distribution of points then the series is mostly stocastic.

For example, for a random series (white noise with mean = 0, standard deviation = 1) the recurrence plot is as shown below (d = 2, h = 1, r = 0.05)

Note that the recurrence plot is diagonally symmetric since the distance of the i'th embedded vector to the j'th embedded vector is the same as the distance of the j'th to the i'th. Also, there is a diagonal line where i = j (the distance between the vectors is 0 and hence always less than the threshold)

 

Long straight parallel lines indicate a periodic series.

For example the following is the recurrence plot for a sinusoid, the distance between the diagonal lines is the period.

Horizontal and vertical white lines or bands are associated with states that don't occur often. Horizontal and vertical black lines indicate states that remain unchanged for periods of time.

 

One of the difficulties is determining an appropriate distance r. In the software developed here and listed below, the user can select a value of r, the software then creates 3 recurrence plots with r/2, r, and 2r. Another approach is to show all radii and to illustrate this with a colour recurrence plot where the range of radii are mapped onto different colours. In the example on the right the same sinusoid is used as in the last example and the radius r is mapped onto a blue to red colour map (blue = small r, red = large r).

 

For a chaotic series on the right the diagram is more complicated. There are brief periodic pieces (short diagonal strips). The time series used below is the logistic equation operating in the chaotic region, namely

x = 4 * xn-1 ( 1 - xn-1 )
 

Significant density changes near the diagonal is an indication of non-stationarity. The following plot is for a linearly swept chirp.

 

The correlation integral is just the "density" of points on a recurrence plot, that is, the number of points divided by the total number of points in the plot.

C(r) = (number of times ||yi - yj|| < r) / N2

Note, the formulation includes the points along the diagonal where i = j. Determining the slope of the correlation integral as a function of log(r) gives the correlation dimension. This is most efficiently computed from a colour recurrence plot where, in general, all the distances have been computed and one simply need to step through a range of radii and count the number of cells with a distance less than the radius.

In the graph on the right the correlation integral is computed over a range of radii for both a white noise source and the logistic equation operating in the chaotic regime. The embedding dimension is 2 and the delay is 1.

 

The graph on the right uses the same chaotic series as above but varies the embedding dimension from 2 to 6 while keeping the delay fixed at 1.

 

In this case the delay is varied from 1 to 3 with a fixed embedding dimension of 2.

 

The usefulness of this result stems from an observation that the fractal dimension of an attractor of a chaotic system exhibits self similarity. It can be shown that the topological dimension of an embedding is the same as the topological dimension of the underlying attractor as long as the embedding dimension is large enough. Unfortunately there is no straightforward way to determine what "large enough" is except by computing the correlation dimension for a number of different embedding dimensions.
Note, the delay used for the embedding is normally taken to be at the first zero of the autocorrelation function.


Further examples

Delta function train  
Sinusoid + noise  
Square wave  
Logistic (double bifurcation)

xn = 3.5 * xn-1 ( 1 - xn-1 )

This example also illustrates another small detail, normally one does not wish to make a recurrence map that is as wide and tall as the number of points in the time series. The time series used in these examples is 2000 points wide while the recurrence maps are 400 pixels wide an tall. This is solved by subsampling but depending on the dimension and delay there can be a border at the edge of the recurrence plot for which the calculations cannot be performed (see the white fringe on the right and top of this recurrence plot).

 
Decaying sinusoid

A change in the thickness of the lines tends to indicate a change in the amplitude of the series.

 

Source code

The following illustrative C source was used to generate the diagrams shown above


References

Understanding Nonlinear Dynamics
Daniel Kaplan and Leon Glass
Springer-Verlag

Recurrence Plots of Dynamical Systems
Eckman, J.P., et. al.
Europhys. Lett. 4, 973 (1 November 1987).

Characterization of Strange Attractors
Grassberger, P., Procaccia, I.
Phys. Rev. Lett. 50, 346 (31 January 1983).