The fundamental starting point of many approaches in nonlinear data analysis is the construction of a phase space portrait of the considered system. The state of a system can be described by its state variables $x^1(t), x^2(t), ... ,x^d(t)$, for example the both state variables temperature and pressure for a thermodynamic system. The d state variables at time t form a vector in a d-dimensional space which is called phase space. The state of a system typically changes in time, and, hence, the vector in the phase space describes a trajectory representing the time evolution, the dynamics, of the system. The shape of the trajectory gives hints about the system; periodic or chaotic systems have characteristic phase space portraits.
The observation of a real process usually does not yield all possible state variables. Either not all state variables are known or not all of them can be measured. However, due to the couplings between the system's components, we can reconstruct a phase space trajectory from a single observation u_i by a time delay embedding (Takens, 1981):
where $m$ is the embedding dimension and $\tau$ is the time delay (index based; the real time delay is $\tau\,\Delta t$). This reconstruction of the phase space is called time delay embedding. The phase space reconstruction is not exactly the same to the original phase space, but its topological properties are preserved, if the embedding dimension is large enough (the embedding dimension has to be larger then twice the phase space dimension, or exactly m > 2 d + 1). And this reconstructed trajectory is sufficient enough for a subsequent analysis.
Now we look at the phase space portrait of an harmonic oscillation, like an undamped pendulum. First we create the position vector y1 and the velocity vector y2
x = 0 : pi/10 : 6 * pi; y1 = sin(x); y2 = cos(x);
The phase space portrait
plot(y1, y2) xlabel('y_1'), ylabel('y_2')
is a circle, revealing the exactly recurrence of each state after one cycle. Using the time delay embedding, we can reconstruct this phase space portrait using only one observation, e.g. the velocity vector, and a delay of 5 (which corresponds to a quarter of the period of our pendulum)
t = 5; plot(y2(1:end - t), y2(1 + t:end)) xlabel('y_2(t)'), ylabel('y_2(t+5)')
As we see, the reconstructed phase space is almost the same as the original phase space.
Next we compare this phase space portrait with such of a typical nonlinear system
x = load('lorenz.dat'); plot(x(:,1), x(:,2)) xlabel('Time (sec)'), ylabel('Temperature')
This dataset contains a time axis in the first column and three state vectors of the Lorenz system in the columns two to four. The Lorenz system is a simplified weather modell describing turbulence with three states (two temperature distributions and velocity).
We plot the phase space portrait with
plot3(x(:,2),x(:,3),x(:,4)) view(-173,42) xlabel('x_1'), ylabel('x_2'), zlabel('x_3')
The first difference to the periodic system above is, that there is not only one circle-shaped trajectory, but a lot of parallel running lines. The trajectory is circling one fixed point, and after some random time, circling another. It is jumping between both circuitries at random instants. The trajectory does not use the same way again, as in the example of the harmonic oscillation, but it recurs to a state visited before very close. Moreover, if we follow two very close segments of the trajectory, we will see that they run into different regions of the phase space after some time. Assuming we have an infinite long trajectory (which is also called attractor), then we can imagine, that this trajectory fills a 2-dimensional plane. But this is not the case; a measurement of the dimension of this trajectory would reveal a fractal dimension between two and three! (We will focus on the dimension estimation of attractors in the next chapter).
These observed properties are typical for chaotic systems: Although small disturbations of such a system cause exponentially divergence of its state, after some time the system will come back to a state that is arbitrary close to a former state and pass through a similar evolution. The infinite long phase space trajectory has a fractal dimension.
The reconstruction of the phase space portrait using only the first state
tau = 15; plot3(x(1:end-2*tau,2),x(1+tau:end-tau,2),x(1+2*tau:end,2)) view(-105,32) xlabel('x_1'), ylabel('x_2'), zlabel('x_3')
reveals a similar phase portrait with the two typical ears. The characteristic properties of chaotic systems can be revealed from this reconstruction as well.
The time delay and embedding dimension has to be chosen with a preceding analysis of the data. The delay can be estimated with the help of the autocovariance or autocorrelation function.
For our example of a periodic oscillation
x = 0 : pi/10 : 6 * pi; y1 = sin(x);
we compute and plot the autocorrelation function
N = length(y1); for i = 1 : 30 r = corrcoef(y1(1 : N - i), y1(1 + i : end)); C(i) = r(1,2); end plot(C) xlabel('Delay'), ylabel('Autocorrelation') grid on
Now we chose such a delay, were the autocorrelation function has its first zero. In our case this is 5, which is the same value which we have already used in our example of phase space reconstruction above. Another (and usually better) measure for chosing the delay is the mutual information (cf. the book "Nonlinear Time Series Analysis" by Kantz and Schreiber).
Using the CRP toolbox, we can compute the mutual information of time delayed time series with the command
Error using ==> evalin Subscript indices must either be real positive integers or logicals.
The mutual information is low where the auto-correlation is low. However, in contrast to the auto-correlation it measures also non-linear dependencies. As delay we use the first minimum of the mutual information or, if such a minimum is not apparent, we use the first delay, where a steep decay of the mutual information changes to a slow decay or a plateau.
The estimation of an appropriate embedding dimension is using the amount of false nearest neighbours. This approach is checking whether the neighbourhoods of all points in the phase space change by increasing the embedding dimension. If the embedding is too low, many neighbours are only projections. Thus, such an embedding dimension would be optimal, where the amount false nearest neighbours (almost) vanishes. Let us look about the amount of false nearest neighbours of the periodic oscillation by using a command from the CRP toolbox.
Here we find a dimension of two as appropriate for the periodic oscillation.
The phase space trajectory or its reconstruction is the base of several measures defined in nonlinear data analysis, like Lyapunov exponents, Renyi entropies or dimensions. It is also the base for the construction of recurrence plots.
F. Takens: Detecting Strange Attractors in Turbulence. In: Dynamical Systems and Turbulence, Lecture Notes in Mathematics, 898, Berlin, 1981, 366-381 H. Kantz, T. Schreiber: Nonlinear Time Series Analysis. University Press, Cambridge, 1997.
1. Phase space trajectory of the Roessler system
1.1. Import the Roessler system from the file roessler.dat and plot the corresponding phase space trajectory.
1.2. Reconstruct the phase space by using the first component and then by the third component. Use different delay parameters. Which is the optimal delay?
2. The file ecg.dat contains data of an ECG measurement.
2.1. Determine an appropriate delay for a phase space reconstruction by using the auto-correlation function.
2.2. Represent the phase space reconstruction for embedding dimension m=3.
2.3. What happens, if the delay is changed?