Getting Started

We start with a time series from the Lorenz equations.

using OrdinaryDiffEq
using Plots, StatsPlots

const σ = 10.0
const ρ = 28.0
const β = 8/3

function lorenz!(du, u, p, t)
    du[1] = σ * (u[2] - u[1])
    du[2] = u[1] * (ρ - u[3]) - u[2]
    du[3] = u[1] * u[2] - β * u[3]
end

prob = ODEProblem(lorenz!, [1.0, 0.0, 0.0], (0.0, 40.0))
sol  = solve(prob)

plot(sol, idxs=(1, 2, 3), xlabel="x", ylabel="y", zlabel="z", legend=false)

Although we usually recommend it, for this example we do not enforce any condition on the sampling density through the integration.

Creating a TrajectorySpace

A TrajectorySpace represents a time series and a cubical cover $Y$. The cycling signature of a segment is the subspace of $H_1(Y)$ induced by including a thickening of the segment into $Y$.

We generate time series data by evaluating the integrated trajectory on an equidistant grid. We then compute the direction (i.e. normalized tangent vector) at each point.

using LinearAlgebra

dt = 0.01
tgrid = 0:dt:40

# time series vector
X = Array(sol(tgrid))

# compute unit tangent vectors
TX = mapslices(X, dims=1) do x
    v = zeros(3)
    lorenz!(v, x, 0, 0)
    return normalize(v)
end

We can construct a TrajectorySpace simply by specifying the time series and two parameters for the size of the boxes used in the cubical cover. More precisely, the cover uses boxes of size boxsize in space and 1/sb_radius in the unit tangent direction (more precisely, the unit tangent vectors are scaled to sb_radius and then covered with boxes with unit side length).

using CyclingSignatures

x_boxsize = 8.0
sb_radius = 1

traj_space = utb_trajectory_space_from_trajectory(X, TX, x_boxsize, sb_radius)

TrajectorySpace{RefinedEquidistantTrajectory{Float64, Matrix{Float64}}, CyclingSignatures.SBCubicalComparisonSpace{CyclingSignatures.CubicalVRCarrier{Int64, SparseArrays.SparseMatrixCSC{Int64, Int64}}, Float64, Int64}, DynamicDistance{Float64}, Float64}(RefinedEquidistantTrajectory 4001 points in 4001 steps., SBCubicalComparisonSpace (boxsize: 8.0, sb_radius: 1, β_1: 2)
  Carrier: CubicalVRCarrier
  Number of points: 166
  Number of edges:  2121
  β_1:              2
, DynamicDistance{Float64}(3, 8.0), 8.0)

Subsegment Experiments

We analyze cycling properties for a collection of randomly sampled segments of different lengths. In the following, we select 100 segments of lengths 10:10:500 from the given dataset.

segment_lengths = collect(10:10:500)
n_runs = 100
exp = RandomSubsegmentExperiment(traj_space, segment_lengths, n_runs, 42)

results = run_experiment(exp)

Result(100 runs for 50 segment lengths)

The resulting collection of cycling signatures can be analyzed to obtain a coarse description of the dynamics on the Lorenz attractor.

evaluation_radius = 3.0

plt_rank = plot_rank_distribution_at_r(results, evaluation_radius)
sig1, plt1 = plot_subspace_frequency_at_r(results, 1, evaluation_radius; n_subspaces=3)
sig2, plt2 = plot_subspace_frequency_at_r(results, 2, evaluation_radius)
pltinc = plot_cycspace_inclusion(sig1, sig2)

combined = plot(plt_rank, plt1, plt2, pltinc; size=(800, 600), layout = (2, 2))

The top-left plot shows the distribution of cycling ranks. Very short segments have rank 0, which is expected because they are too short to wrap around anything. After a certain time span, most segments have nontrivial cycling rank, indicating that cycling is typical in this time series.

The subspace-frequency plots summarize which cycling spaces occur most often. In the Lorenz system, rank-1 signatures correspond to the dominant cycling motions around the two wings and transitions between them. Longer datasets and denser segment-length grids reveal more fine-grained recurrent structure.