Recurrence Plots

Recurrence Matrices

A Recurrence plot (which refers to the plot of a matrix) is a way to quantify recurrences that occur in a trajectory. A recurrence happens when a trajectory visits the same neighborhood on the phase space that it was at some previous time.

The central structure used in these recurrences is the (cross-) recurrence matrix:

where $d(x[i], y[j])$ stands for the distance between trajectory $x$ at point $i$ and trajectory $y$ at point $j$. Both $x, y$ can be single timeseries, full trajectories or embedded timeseries (which are also trajectories).

If $x\equiv y$ then $R$ is called recurrence matrix, otherwise it is called cross-recurrence matrix. There is also the joint-recurrence variant, see below. With RecurrenceAnalysis you can use the following functions to access these matrices

RecurrenceMatrix(x, ε; kwargs...)

CrossRecurrenceMatrix(x, y, ε; kwargs...)

JointRecurrenceMatrix(x, y, ε; kwargs...)

Simple Recurrence Plots

The recurrence matrices are internally stored as sparse matrices with boolean values. Typically in the literature one does not "see" the matrices themselves but instead a plot of them (hence "Recurrence Plots"). By default, when a Recurrence Matrix is created we "show" a mini plot of it which is a text-based scatterplot.

Here is an example recurrence plot/matrix of a full trajectory of the Roessler system:

using DynamicalSystems
ro = Systems.roessler(ones(3), a=0.15, b=0.20, c=10.0)
N = 2000; dt = 0.05
tr = trajectory(ro, N*dt; dt = dt, Ttr = 10.0)

R = RecurrenceMatrix(tr, 5.0; metric = "euclidean")
recurrenceplot(R; ascii = true)

      RecurrenceMatrix of size (2001, 2001) with 383597 entries
    +------------------------------------------------------------+ 
    |         ..:'.::'       .::'.:'..:'       ..:'.::    ..  .:'|
    |   ..  .::'.:''   .:' .:'..''.::'   .   .::'.:'    .:'..:'  |
    | .:':.:':::' .  .:'.::'.::'.:'    .:'.::'.::'..  ::'.::'.:' |
    |:'.::'.::'      :.:::.:'    ..    '.::'.::'      :.:::.:'   |
    | ''' ''' ' ..  .:'''''   ..:'..  .:'''''  ' .. ..:' '''   ..|
    |         .:'.::'       .::'.:'.::'       .::'.::'    .  .:''|
    |      .:' .:'..:    ..:'..:'.:''.:'    .:'..:'..:    ..:'.::|
    |   ..:'.::'.::'   .::'.::' :'.::'   .::'.::'.:'' ' .::'.:'' |
    |'.::'.::'.:'    .:::.:'..  .:''   .:::.:::.:'    .:':::'    |
    |:::::':::'      '.::'.:''  '      ':::'.::'      '.::'..    |
    |'::' ::'   ..  .::'::'    .. .   .::'.::'   ..  .::'::'    .|
    |         .:'.::'       .::'.:' .:'       ..:'.::'       .::'|
    |      .::'.::'.:     .:''.:'.::'..     .::'.:''.: .   .:':.:|
    |     :'' :'.::'.:  ::' ::'.:':.:'..   :' ::'.::'.:' ::'.::'.|
    |     .   .:'..:'    .:'.::'.:''.:'    .  ..:'.::'   ..:'.::'|
    |  ..:'.::'.::'   .::'.::'..'.::'   ..:'.::'.:''   .::'.:''  |
    |.::'.:::.:'     :::.:'.::' ::'    :::.:::.:'     :':::'.::  |
    |::::':::'       .::'.::'   .      .::'.::'       .::'.:''   |
    | ''  ''  '.:' .:'' ''    .:'..: .:'' '''   .:' .:'' ''    .:|
    |        .:'..:'.      .::'.:''.:'..      .:'.::'.    . .::'.|
    |     .::'.::'.:'    .:''.:'.::'.:'    .::'.:':.:'.   .:'..:'|
    |.  .:':.:'..:'  '.::'.::'  ''.:'   ..:':::'..''  '.::'.::'  |
    |:::'.::'.:''    :::.::'..  ::'    ::'.::'.:'     :::.:''    |
    |::::::::'       .::'.::'          .::::::'       .::'.::    |
    |::::::'..       ::.::'            :::::'.        ::.:'      |
    +------------------------------------------------------------+ 

typeof(R)

RecurrenceMatrix

summary(R)

"RecurrenceMatrix of size (2001, 2001) with 383597 entries"

The above simple plotting functionality is possible through the package UnicodePlots. The following function creates the plot:

recurrenceplot([io,] R; minh = 25, maxh = 0.5, ascii, kwargs...) -> u

Here is the same plot but using Unicode Braille characters

recurrenceplot(R; ascii = false)

      RecurrenceMatrix of size (2001, 2001) with 383597 entries
    ┌────────────────────────────────────────────────────────────┐ 
    │⠀⠀⠀⠀⠀⠀⠀⠀⠀⢀⣠⠞⠃⢀⣴⠟⠁⠀⠀⠀⠀⠀⠀⠀⢀⣴⠟⠁⣠⡶⠋⢀⣤⠞⠁⠀⠀⠀⠀⠀⠀⠀⢀⣤⠞⠁⣀⡴⠟⠀⠀⠀⠀⢀⣤⠀⠀⣀⡴⠛│
    │⠀⠀⠀⢀⡀⠀⠀⣀⣴⠟⠁⣠⡶⠋⠁⠀⠀⠀⣀⡴⠋⠀⣠⡾⠋⢀⣠⠘⠉⣀⡴⠟⠁⠀⠀⠀⢀⠀⠀⠀⣠⡶⠟⠁⣠⡾⠋⠀⠀⠀⠀⣠⡶⠋⢀⣠⡾⠋⠀⠀│
    │⠀⣠⡶⠛⢁⣠⡾⠋⢁⣴⠾⠋⠀⡀⠀⠀⢠⡾⠋⢀⣴⠟⠋⣀⡴⠟⠁⢠⡾⠋⠀⠀⠀⠀⢠⡾⠋⢀⣴⡾⠋⢀⣴⠞⠉⢀⡀⠀⠀⣴⡾⠋⣀⣴⠟⠉⣀⡴⠋⠀│
    │⡿⠋⣠⣴⡿⠋⣠⣶⠟⠁⠀⠀⠀⠀⠀⠀⢈⣠⣾⠟⢁⣠⡾⠋⠀⠀⠀⠀⢀⡀⠀⠀⠀⠀⠈⣠⣶⠟⠋⣠⣾⠟⠁⠀⠀⠀⠀⠀⠀⢁⣤⡾⠟⢁⣠⡾⠋⠀⠀⠀│
    │⠀⠚⠛⠁⠀⠚⠛⠁⠀⠒⠀⣠⡄⠀⠀⣠⠞⠛⠁⠐⠛⠋⠀⠀⠀⢀⣠⠞⠋⢀⡄⠀⠀⣠⡼⠛⠁⠐⠛⠋⠀⠀⠒⠀⣠⡄⠀⢀⣤⠟⠋⠀⠐⠛⠉⠀⠀⠀⢀⣤│
    │⠀⠀⠀⠀⠀⠀⠀⠀⠀⣤⠞⠉⢀⣴⠟⠁⠀⠀⠀⠀⠀⠀⠀⢀⣴⠟⠁⣠⡾⠋⢀⣴⠞⠉⠀⠀⠀⠀⠀⠀⠀⢀⣴⠞⠁⣀⡴⠟⠁⠀⠀⠀⠀⡀⠀⠀⣠⡴⠛⠁│
    │⠀⠀⠀⠀⠀⠀⣀⡴⠂⠀⣠⡾⠋⢀⣠⠞⠀⠀⠀⠀⢀⣠⡾⠋⢀⣤⠞⠁⣠⡴⠛⠁⣠⡶⠃⠀⠀⠀⠀⣠⡴⠋⢀⣠⡾⠋⢀⣤⠞⠀⠀⠀⠀⢀⣤⠾⠋⢀⣴⠞│
    │⠀⠀⠀⢀⣤⡾⠋⢀⣴⠟⠋⢀⡴⠟⠁⠀⠀⠀⣀⣴⠟⠁⣠⡴⠟⠁⠀⠸⠋⢀⣴⠞⠁⠀⠀⠀⢀⣴⠾⠋⣀⣴⠟⠁⣀⡴⠋⠁⠀⠒⠀⣠⣴⠟⠁⣠⡶⠋⠁⠀│
    │⠃⣠⣶⠟⠉⣠⣾⠟⠁⣠⡾⠋⠀⠀⠀⠀⣠⣾⠟⢁⣤⡾⠋⢀⡀⠀⠀⢠⡴⠋⠁⠀⠀⠀⢠⣶⠟⢁⣠⡾⠟⢁⣠⠾⠋⠀⠀⠀⠀⣤⡾⠛⢁⣴⡾⠋⠀⠀⠀⠀│
    │⣿⠟⣁⣴⡿⠋⣁⣴⠟⠋⠀⠀⠀⠀⠀⠀⠉⣠⣶⠿⠋⣠⡴⠛⠁⠀⠀⠈⠀⠀⠀⠀⠀⠀⠈⣡⣴⡿⠋⣠⣴⠟⠁⠀⠀⠀⠀⠀⠀⠉⣠⣾⠟⠁⣠⡄⠀⠀⠀⠀│
    │⠁⠾⠿⠋⠀⠾⠟⠁⠀⠀⠀⢀⡀⠀⠀⣠⡾⠟⠁⠰⠾⠋⠀⠀⠀⠀⣠⡄⠀⢀⠀⠀⠀⣀⡼⠿⠋⠠⠾⠟⠁⠀⠀⠀⣀⡀⠀⠀⣠⠿⠟⠁⠰⠞⠋⠀⠀⠀⠀⣠│
    │⠀⠀⠀⠀⠀⠀⠀⠀⠀⣠⡾⠋⢀⣴⠞⠉⠀⠀⠀⠀⠀⠀⠀⢀⣴⠞⠁⢠⡴⠋⠀⣠⡾⠋⠀⠀⠀⠀⠀⠀⠀⢀⣠⠞⠋⢀⣴⠟⠁⠀⠀⠀⠀⠀⠀⠀⢀⣴⠟⠁│
    │⠀⠀⠀⠀⠀⠀⢀⣴⠞⠁⣠⡴⠟⠁⣠⠶⠀⠀⠀⠀⠀⣠⡴⠋⠁⣠⡾⠋⢀⣴⠟⠁⣀⡤⠀⠀⠀⠀⠀⣀⣴⠟⠁⣠⡶⠋⠁⣠⠖⠀⢀⠀⠀⠀⣠⡾⠋⢁⣠⠞│
    │⠀⠀⠀⠀⠀⠶⠋⠁⠀⠾⠋⢀⣴⠞⠁⣠⡴⠀⠀⠰⠾⠋⠀⠰⠟⠁⣠⡼⠛⢁⣤⠞⠋⢀⡄⠀⠀⠀⠾⠋⠀⠰⠾⠋⣀⡴⠟⠁⣠⠾⠋⠀⠰⠟⠋⢀⡰⠟⠁⣠│
    │⠀⠀⠀⠀⠀⢀⠀⠀⠀⣠⡾⠋⢀⣤⠞⠋⠀⠀⠀⠀⣠⡶⠃⢀⣴⠞⠋⢠⡴⠛⠁⣠⡾⠋⠀⠀⠀⠀⣀⠀⠀⢀⣠⡾⠋⢀⣴⠞⠁⠀⠀⠀⢀⣠⠞⠋⢀⣴⠞⠁│
    │⠀⠀⢀⣠⡾⠋⢀⣴⠟⠋⢀⣴⠟⠁⠀⠀⠀⣀⣴⠟⠁⣠⣴⠟⠁⣠⠄⠈⢀⣴⠞⠉⠀⠀⠀⢀⣤⠞⠉⣀⣴⠟⠁⣠⡴⠛⠁⠀⠀⠀⣀⣴⠟⠁⣠⡶⠛⠁⠀⠀│
    │⣠⣴⠟⠁⣠⣾⠟⢁⣠⡾⠋⠀⠀⠀⠀⠀⠸⠟⢁⣤⡾⠋⢀⣴⠞⠁⠀⠰⠟⠁⠀⠀⠀⠀⠰⠟⢁⣠⡾⠟⢁⣤⡾⠋⠀⠀⠀⠀⠀⠾⠋⢁⣴⡾⠋⢀⣴⠆⠀⠀│
    │⠟⣁⣴⡿⠛⣁⣴⠿⠋⠀⠀⠀⠀⠀⠀⠀⢠⣶⡿⠋⣠⣶⠟⠁⠀⠀⠀⢀⠀⠀⠀⠀⠀⠀⢠⣴⡿⠋⣠⣴⠟⠉⠀⠀⠀⠀⠀⠀⠀⣠⣶⠟⠉⣠⡾⠋⠁⠀⠀⠀│
    │⠀⠉⠉⠀⠀⠉⠁⠀⠀⠁⢀⡴⠂⠀⣠⡶⠋⠁⠀⠈⠉⠀⠀⠀⠀⣠⡶⠋⢀⣤⠆⠀⣀⡴⠋⠉⠀⠈⠉⠁⠀⠀⠀⣀⡴⠂⠀⣠⡶⠋⠁⠀⠈⠉⠀⠀⠀⠀⣠⠶│
    │⠀⠀⠀⠀⠀⠀⠀⠀⣠⠾⠋⢀⣤⠞⠋⢀⠀⠀⠀⠀⠀⠀⢀⣴⠞⠉⣀⡰⠋⠁⣠⡾⠋⢀⡀⠀⠀⠀⠀⠀⠀⣠⠾⠋⢀⣴⠟⠁⢀⠀⠀⠀⠀⠄⠀⢀⣴⠟⠁⣀│
    │⠀⠀⠀⠀⠀⢀⣴⠞⠁⣠⡴⠟⠁⣠⡶⠋⠀⠀⠀⠀⣠⡶⠛⠁⣠⡾⠋⢀⣴⠟⠁⣀⡴⠋⠀⠀⠀⠀⣀⣴⠟⠁⣠⡶⠋⢁⣠⠾⠋⢀⠀⠀⠀⣠⡾⠋⢀⣠⠾⠋│
    │⡤⠀⠀⣠⡾⠛⢁⣤⡾⠋⢀⣤⠞⠁⠀⠀⠈⢀⣴⡾⠋⢀⣴⠟⠉⠀⠀⠈⠁⣠⠾⠋⠀⠀⠀⢀⣠⡾⠋⢁⣴⡾⠋⢀⡤⠈⠁⠀⠀⠉⢀⣴⠿⠋⣀⣴⠟⠁⠀⠀│
    │⣁⣴⡿⠋⣠⣴⠟⠋⣠⡴⠛⠁⠀⠀⠀⠀⣼⠟⢉⣠⡾⠟⠁⣠⠄⠀⠀⢰⠟⠁⠀⠀⠀⠀⢰⡿⠋⣠⣶⠟⠁⣠⡶⠋⠀⠀⠀⠀⠀⣾⠟⢁⣠⡾⠋⠁⠀⠀⠀⠀│
    │⡿⢋⣴⣾⠟⢁⣴⡾⠋⠀⠀⠀⠀⠀⠀⠀⢠⣴⡿⠋⣀⣴⠟⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⢀⣴⣿⠟⣁⣴⡿⠋⠀⠀⠀⠀⠀⠀⠀⣠⣴⡿⠋⣀⣴⠆⠀⠀⠀⠀│
    │⣾⣿⢟⣡⣾⡿⠋⣀⡀⠀⠀⠀⠀⠀⠀⠀⠻⢋⣤⣾⠟⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠸⢋⣥⣾⠟⠉⣀⠀⠀⠀⠀⠀⠀⠀⠀⠟⢋⣤⡾⠋⠀⠀⠀⠀⠀⠀│
    └────────────────────────────────────────────────────────────┘ 

As you can see, the Unicode based plotting doesn't display nicely everywhere. It does display perfectly in e.g. Juno, which is where it is the default printing type. Here is how it looks like in a dark background:

Advanced Recurrence Plots

A text-based plot is cool, fast and simple. But often one needs the full resolution offered by the data of a recurrence matrix.

There are two more ways to plot a recurrence matrix using RecurrenceAnalysis:

coordinates(R) -> xs, ys

grayscale(R [, bwcode]; width::Int, height::Int, exactsize=false)

For example, here is the representation of the above R from the Roessler system using both plotting approaches:

using PyPlot
figure(figsize = (10,5))

ax = subplot(121)
xs, ys = coordinates(R)
scatter(xs, ys, color = "k", s = 1)
xlim(1, size(R)[1]); ylim(1, size(R)[2]);
ax.set_aspect("equal")

subplot(122)
Rg = grayscale(R)
imshow(Rg, cmap = "binary_r", extent = (1, size(R)[1], 1, size(R)[2]))

and here is exactly the same process, but using the embedded trajectory instead

y = tr[:, 2]
τ = estimate_delay(y, "mi_min")
m = reconstruct(y, 2, τ)
R = RecurrenceMatrix(m, 5.0; metric = "euclidean")

figure(figsize = (5,5))

xs, ys = coordinates(R)
scatter(xs, ys, color = "k", s = 1)
xlim(1, size(R)[1]); ylim(1, size(R)[2]);

which justifies why recurrence plots are so fitting to be used in embedded timeseries.

Example

In the following we will plot recurrence plots of the Lorenz system for a periodic and chaotic regime (using scatter plot).

lor = Systems.lorenz()
figure(figsize = (10,10))

for (i, ρ) in enumerate((69.75, 28.0))
    set_parameter!(lor, 2, ρ)
    t, dt = 20.0, 0.01
    tr = trajectory(lor, t; dt = dt, Ttr = 2000.0)
    tvec = 0:dt:t

    subplot(2,2, i)
    plot(tr[:, 1], tr[:, 3], color = "C$(i+1)", label = "X vs Z")
    title("ρ = $ρ, " * (i != 1 ? "not periodic" : "periodic")); legend()

    ε = i == 1 ? 5.0 : 3.0
    R = RecurrenceMatrix(tr, ε)

    subplot(2,2,i+2)
    x, y = coordinates(R)
    scatter(tvec[x], tvec[y], s = 1, alpha = 0.2, color = "C$(i+1)")
    xlim(0, t); ylim(0, t); gca().set_aspect("equal")
    xlabel("t"); i == 1 && ylabel("t");
end
PyPlot.tight_layout()

┌ Warning: Assignment to `dt` in soft scope is ambiguous because a global variable by the same name exists: `dt` will be treated as a new local. Disambiguate by using `local dt` to suppress this warning or `global dt` to assign to the existing global variable.
└ @ none:4
┌ Warning: Assignment to `tr` in soft scope is ambiguous because a global variable by the same name exists: `tr` will be treated as a new local. Disambiguate by using `local tr` to suppress this warning or `global tr` to assign to the existing global variable.
└ @ none:5
┌ Warning: Assignment to `R` in soft scope is ambiguous because a global variable by the same name exists: `R` will be treated as a new local. Disambiguate by using `local R` to suppress this warning or `global R` to assign to the existing global variable.
└ @ none:13
┌ Warning: Assignment to `y` in soft scope is ambiguous because a global variable by the same name exists: `y` will be treated as a new local. Disambiguate by using `local y` to suppress this warning or `global y` to assign to the existing global variable.
└ @ none:16

On the left we see long (infinite) diagonals repeated over and over for different times. This is the case for periodic systems as they visit exactly the same area on the phase space again and again. The distance between the offset diagonals also coincides with the periodicity of the system, which is around t ≈ 4.

On the right we see a structure typical of chaotic motion on a strange attractor such as the one of the Lorenz system: the orbit visits neighborhoods of previous points but then quickly diverges again. This results in many small diagonal lines.

Distances

The distance function used in RecurrenceMatrix and co. can be specified either as a string or as any Metric instance from Distances. In addition, the following function returns a matrix with the cross-distances across all points in one or two trajectories:

distancematrix(x [, y = x], metric = "euclidean")