CRP Toolbox

Don't hesitate to ask us questions not included here by using the Recurrence Plot Forum. Please understand that we cannot answer every eMail regarding the CRP toolbox.

I have tried the demo version of the CRP toolbox. Is the full version of the toolbox freely available?

Yes, of course. The toolbox is freely available. However, you have to request an access code, which is necessary to download the toolbox. Moreover, you have to cite this toolbox in your reports, articles and papers, if you use results computed with the toolbox.

top

I have downloaded the install.m file. But how to install the toolbox?

The file install.m is the repository of the toolbox and contains everything of the toolbox. The toolbox will be installed when you execute the install.m file from within MATLAB.

To do so, please open MATLAB, go to the command window, then change to the directory in which the downloaded file install.m is located, and finally type

install

in the MATLAB command window.

This will then automatically unpack the toolbox from the install.m file and copy all files to the standard MATLAB toolbox folder.

top

Why I get a checksum error during installation?

The installation script checks itself for integrity. Only this guarantees a proper installation and usage of the programmes in the toolbox.

Checksum errors may occur if you have received the installation script by eMail, and your eMail client has processed the attachement e.g. with an anti-virus scanner. Please ensure, that eMail attachements will not be changed by anti-virus scanners.

Another possebility is that you have used a proxy server during downloading the toolbox from the web and the proxy server has modified the data in any way. You should try to download the toolbox without a proxy connection.

top

Why I get an error during the installation that the access to a directory is denied ot that a directory cannot be created (MS Windows)?

Under a MS Windows operating system, MATLAB installs toolboxes in the MATLAB path. Usually you will need administrator rights to install files in this path.

A future version of the installation routine will fix this issue.

top

I cannot install but get an errorcode 95.01 (MS Windows).

If you get an errorcode 95.01 and some error related to the userpath function, please check whether the variable userpath is set and whether you have a MATLAB directory in your Documents path (something like C:\Documents and Settings\YOURNAME\Documents\MATLAB). Older MATLAB versions (reported e.g. for MATLAB 2006b) may not be able to reset the userpath variable. Therefore, just create manually in your Documents path a new directory called MATLAB. The installation should then work. If it will not work either, please let us know.

top

I cannot install but get an errorcode 95.02 (textscan error).

If you get an errorcode 95.02 and some error related to the textscan function, please check whether the variable userpath is set. It should be a standard MATLAB variable pointing to the default path where MATLAB looks for user toolboxes. If the variable is not set, then please refer to https://www.mathworks.com/help/matlab/ref/userpath.html and set (or reset) it. If it will not work either, please let us know.

top

I cannot install but get an errorcode 95.02 (invalid directory or directory does not exist).

If you get an errorcode 95.02 and an error invalid directory or directory does not exist, probably the standard MATLAB user path does not exist in your system and MATLAB is not able to create it. The user path is contained in MATLAB variable userpath. Create a folder that corresponds to this path. It might be ~/Documents/MATLAB or ~/matlab.

top

Why I get an error during installation the toolbox?

We have listen of rather rare problems during installation the toolbox. (e.g. a syntax error may occur regarding a missing rehash function in MATLAB, or the version command in MATLAB obtains strange (and non-standard) version information). Usually these errors are due to an out-dated MATLAB version (before release 11 or version 5.3, resp.). You should consider to update MATLAB to a newer version.

top

Starting MATLAB, I get a strange error about failed XML validation. What happened?

The toolbox contains an XML file which gives some information to the MATLAB environment. However, Mathworks has changed the style of this XML file with each new release of MATLAB. Therefore, I gave up to fix this issue for some time. With release 26.2 of the CRP Toolbox, this issue was fixed. If you are bored with these stupid error message, either just remove all XML files from the toolbox folder (if you don't know the location of the toolbox, just type which CRPtool on the MATLAB commandline) or install a newer version of the CRP Toolbox. Nevertheless, the toolbox is working fine and without any problem despite the XML bug.

top

When starting any function from the toolbox, I get the error "Undefined function 'narginchk'."

Unfortunately, Mathworks decided to break upwards and downwards compatibility. The previous version for checking the input arguments of a function will be removed in a future version of MATLAB. However, the new command for input argument validation is only available since MATLAB 2011b. To solve this problem, you can either update to a newer version of MATLAB or install a custom made narginchk function in your MATLAB folder, e.g., by using the freely available code at https://gist.github.com/hagenw/5642886. It is suggested to copy this file into your toolbox folder which can be located by the command which crp (the toolbox path is everything before "CRPtool").

top

I'm using the rp_plugin on a PC with MS Windows. Why the plugin is not working for computation of recurrence plots?

There is still a problem with the plugin for the MS Windows. If you are using the plugin to create a recurrence plot with, e.g., crp, it is not working. We are not using windows PCs, therefore, we have not a big chance to find out the reason. Sorry.

Nevertheless, the RQA can be computed by the plugin.

Other operating systems are not concerned.

top

Where can I find a documentation on how to properly use the CRP toolbox/ recurrence analysis?

A printable documentation for the proper use of the CRP toolbox is available under the subsection printable reference manual.

A brief introduction and a forum about recurrence plots, cross and joint recurrence plots as well as their quantification can be found at www.recurrence-plot.tk.

top

How do I cite the toolbox in my paper?

The use of the CRP Toolbox should be cited as

• N. Marwan: Cross Recurrence Plot Toolbox for MATLAB®, Ver. 5.24 (R34), https://tocsy.pik-potsdam.de/CRPtoolbox/, accessed 2022-06-29.

Moreover, at least the following paper should also be cited

• N. Marwan, M. C. Romano, M. Thiel, J. Kurths: Recurrence Plots for the Analysis of Complex Systems, Physics Reports, 438(5–6), 237–329 (2007). DOI:10.1016/j.physrep.2006.11.001

Please check the User Agreement for the conditions of using this software and which publications should be cited in your work.

top

How long the data series can be?

This depends on your computer power, especially on the memory, and the kind of computation needed. For higher speed in output the whole matrix of the (cross/ joint) recurrence plot is in the work space – which is limited by the available computer memory (unfortunately, MATLAB is rather memory hungry). On current computers, the maximal data length can be up to 5,000 (could be longer, test it on your system) for crqa of the entire data (window length equals data length). However, using rqa the length of used data can be much larger using small windows.

For recurrence plots, using crp the data length can be up to 5,000 (could be longer, test it on your system). For longer data sets, crp_big is more appropriate. Moreover, read How can I compute a recurrence plot for really long data series? if you really need a (cross/ joint) recurrence plot for long data series.

If your system is supported, you can use a plugin for the CRP Toolbox, which allows to compute (cross/ joint) recurrence plots and their quantification for longer data series. The used data can be more than 10 times longer than without using the plugin!

Alternatively, you can use a commandline programme for computing recurrence plots and RQA (without MATLAB).

top

How can I compute a recurrence plot for really long data series?

For computing (cross/ joint) recurrence plots of long data series, use a similar script as in the following. The data length is finally limited by the used platform performance. The examples also illustrate the capability of using the programms in a script. The first example uses sparse matrices (when calculating a recurrence plot of only one time series, replace x2 by x1):

%% some parameter settings and create example time series
m = 3; t = 20; e = .5; w = 300;

x1 = sin((1:5000)/40)'; x2 = sin((1:7000)/80)';

Y = spalloc(length(x2)-(m-1)*t,length(x1)-(m-1)*t,1);
k = 0; h1 = waitbar(0,'Compute sub CRPs - Please be patient.');
Nx = length(x1)-(m-1)*t; Ny = length(x2)-(m-1)*t;
ax = ceil(Nx/w); ay = ceil(Ny/w);
Nx2 = floor(Nx/ax); Ny2 = floor(Ny/ay);

%% compute single CRPs and fill the sparse matrix
for i = 1:Nx2:Nx-Nx2;
for j = 1:Ny2:Ny-Ny2, k = k+1; waitbar(k/(Nx*Ny/(Nx2*Ny2)))

X2 = crp(x1(i:i+Nx2+(m-1)*t), x2(j:j+Ny2+(m-1)*t),m,t,e,...
'nonorm','max','silent');
X = sparse(double(X2));
Y(j:j+Ny2-1, i:i+Nx2-1) = X(1:Ny2,1:Nx2);

end
end
close(h1)

spy(Y)


The second example writes single (cross/ joint) recurrence plots to the hard disk (when calculating a recurrence plot of only one time series, replace x2 by x1):

%% some parameter settings and create example time series
m = 3; t = 20; e = .5; w = 300;

x1 = sin((1:5000)/40)'; x2 = sin((1:7000)/80)';
Nx = length(x1); Ny = length(x2);

%% compute single CRPs and write them to the hard disk
b1 = zeros((m-1)*t+ceil(length(x1)/w)*w,1);
b1(1:length(x1)) = x1;
b2 = zeros((m-1)*t+ceil(length(x2)/w)*w,1);
b2(1:length(x2)) = x2;
h = waitbar(0,'Compute sub CRPs - Please be patient.')

for i = 1:w:length(b1)-w-1, waitbar(i/((length(b1)-w-1)))
for j = 1:w:length(b2)-w-1
X = crp(b1(i:i+w+(m-1)*t-1),b2(j:j+w+(m-1)*t-1),m,t,e,...
'max','silent','nonorm');
i2 = num2str((i+w-1)/w);j2=num2str((j+w-1)/w);

filename = ['CRP_',i2,'_',j2,'.tif'];
imwrite(X,filename,'tif')
end
end, close(h)

%% read single CRPs and unify them
xmax = (i+w-1)/w; ymax=(j+w-1)/w;
Y = zeros(length(b1),length(b2));

for i = 1:xmax,waitbar(i/xmax)
for j = 1:ymax
i2 = num2str(i);j2=num2str(j);
filename = ['CRP_',i2,'_',j2,'.tif'];

Y(i*w-(w-1):i*w,j*w-(w-1):j*w) = X';
end
end, close(h)

Y(Nx+1:end,:) = [];Y(:,Ny+1:end) = [];
spy(double(Y))


Alternatively, you can use a plugin for the CRP Toolbox, or a commandline programme for computing recurrence plots (without MATLAB).

top

Will the data be scaled by the programmes?

The data will be scaled by default, i.e. the data will be normalized to a mean of zero and a standard-deviation of one. This is a normal procedure in data analysis. This may also be helpful for applying CRPs to two data series with different magnitudes. However, you can use the programmes without normalization the data. Just call the programmes with the argument 'nonormalize', e.g.

crp(x,y,3,4,.1,'nonormalize')
or shorter
crp(x,y,3,4,.1,'non')


top

Does it be possible to use original phase space vectors instead of embedding?

Yes, of course. The commands crp2 and crqa support original phase space vectors. E.g., creation of RP from the phase space vectors of a harmonic pendulum:

x = [sin(0.1 * [1:100])', cos(0.1 * [1:100])'];
crp2(x)


top

How to chose an appropriate threshold value?

There are different possibilities or requirements for the threshold.

First, you should consider the articles

• M. Thiel, M. C. Romano, J. Kurths, R. Meucci, E. Allaria, F. T. Arecchi: Influence of observational noise on the recurrence quantification analysis, Physica D, 171(3), 138–152 (2002). DOI:10.1016/S0167-2789(02)00586-9
• S. Schinkel, O. Dimigen, N. Marwan: Selection of recurrence threshold for signal detection, European Physical Journal – Special Topics, 164(1), 45–53 (2008). DOI:10.1140/epjst/e2008-00833-5
• K. H. Kraemer, R. V. Donner, J. Heitzig, N. Marwan: Recurrence threshold selection for obtaining robust recurrence characteristics in different embedding dimensions, Chaos, 28(8), 085720 (2018). DOI:10.1063/1.5024914

A too large threshold makes usely no real sense. Mostly it is suggested to choose the threshold in such a way, that its value corresponds to 10% of the maximum or mean phase space diameter (using the CRP toolbox you can find these diameters with the command pss). Another way is to choose the threshold in such a way, that the recurrence rate is 10%. Further possibilities are to compute several RQA measures for an increasing threshold value and to look for a region where the RQA measures change slowly.

top

What is the meaning of embedding dimension, time delay, and how to chose them?

The state of a system can be described by its $$d$$ state variables $$x_1(t), x_2(t), \ldots, x_d(t)$$ and can be written as a $$d$$-dimensional vector $$\vec{x}(t)$$. However, 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. Most often only one observation (measurement) $$u(t) = u_k$$ (with $$i = \Delta t$$) is available. Following Takens' embedding theorem (1981) we can reconstruct the phase space from a single time series $$u_k$$ by using an embedding dimension $$m$$ and a time delay $$\tau$$

$$\vec{x}(t) = \vec{x}_i = ( u_i, u_{i+\tau}, \ldots, u_{i+(m-1)\tau} ), \qquad t = i \Delta t,$$

whereas $$\vec{x}(t)$$ is the vector of reconstructed states in the phase-space at the time $$t$$. The choice of $$m$$ and $$\tau$$ should base on methods for detecting the optimal values of these two parameters like method of false nearest neighbours, fnn (for $$m$$) and mutual information, mi (for $$\tau$$), which ensures the entire covering of all free parameters and avoiding autocorrelated effects (e.g. Kantz and Schreiber, 1997).

• N. Marwan: Encounters With Neighbours – Current Developments Of Concepts Based On Recurrence Plots And Their Applications, Ph.D. Thesis, University of Potsdam, ISBN 3-00-012347-4, urn:nbn:de:kobv:517-0000856, 2003

top

What is the meaning of "vector switching" in the control panel of crp, crp_big, crp2 and jrp?

Vector switching means, that single components of the phase space vector will get a negative sign. For instance, switching the 2nd component of a state $$\vec{x}(t)$$ of a 3-dimensional system:

$$\left( \begin{array}{rrr} x_1(t)\\x_2(t)\\ x_3(t)\\ \end{array} \right) \rightarrow \left( \begin{array}{rrr} x_1(t)\\ -x_2(t)\\ x_3(t)\\ \end{array} \right)$$.

top

What is the meaning of the different norms and recurrence criteria (maximum, Euclidean, normalized, fixed, order patterns etc.)?

Different norms and recurrence criteria (neighbourhood criteria) can be used for the definition of the recurrence of a state.

 Maximum Norm ($$L_\infty$$-norm) the distance between two phase space vectors is the maximal distance between their components $$\max\bigl(|x_1(t_1) - x_1(t_2)|, \ldots, |x_m(t_1) - x_m(t_2)|\bigr)$$ Euclidean Norm ($$L_2$$-norm) the distance between two phase space vectors is the smallest distance between the both vectors $$\sqrt {\sum_i (x_i(t_1) - x_i(t_2))^2}$$ Minimum Norm ($$L_1$$-norm, Manhattan or Taxicab norm) the distance between two phase space vectors is the sum of the distances of all components $$\sum_i |x_i(t_1) - x_i(t_2)|$$ Normalized Norm all phase space vectors are normalised to a length of one $$\vec{x}'(t) = \frac{\vec{x}(t)}{\|\vec{x}(t)\|}$$ then the Euclidean distance is used Fixed amount ofnearest neighbours The number of neighbours in the neighbourhood is constant, i. e. the number of recurrence points in one column of the RP is constant. Such an RP is not symmetric. This is the original definition of an RP by Eckmann et al. (1987). J.-P. Eckmann, S. Oliffson Kamphorst, D. Ruelle: Recurrence Plots of Dynamical Systems, Europhysics Letters, 4(9), 973–977 (1987). DOI:10.1209/0295-5075/4/9/004 Interdependent neighbours Experimental recurrence criterion for CRPs (does not make sense for auto RPs) and is similar to a JRP. A recurrence point is defined, if the indices of the neighbours of the first trajectory coincides with the indices of the neighbours of the second trajectory. The neighbours of the second trajectory are found by a dynamic threshold, given by the radius of the current neighbourhood of the first trajectory. This was inspired by A. Schmitz: Measuring statistical dependence and coupling of subsystems, Physical Review E, 62(5), 7508-7511 (2000). Order matrix A recurrence is defined by the condition $$x_i \ge x_j$$ C. Bandt: Ordinal time series analysis, Ecological Modelling, 182(3–4), 229–238p. (2005). DOI:10.1016/j.ecolmodel.2004.04.003 Order patterns A recurrence is defined by the recurrence of order patterns πi, which present the dynamics of the time series by a symbolisation of local rank orders $$x_i \rightarrow \pi_k \quad$$ where   $$\pi_k \in \mathbb{R}^m$$ A. Groth: Visualization of coupling in time series by order recurrence plots, Physical Review E, 72(4), 046220p. (2005). DOI:10.1103/PhysRevE.72.046220

top

When I use the recurrence criteria of fixed amount of nearest neighbours or order patterns, does it matter to normalise or not the data?

For recurrence plots basing on fixed amount of nearest neighbours (FAN) or order patterns, it doesn't matter to normalise the data, because normalisation doesn't change the neighbourhood (FAN) or the local rank order (order patterns). For cross recurrence plots (CRPs), normalisation affects the result for FAN based CRP, because the two different trajectories in the phase space will be stretched or compressed and, hence, change their intersecting closeness. For cross order patterns recurrence plots it again doesn't matter.

top

Why the RQA results differ from such gained with the RQA software of Charles Webber?

At first we must say, that the CRP toolbox was independently developed from the work of Charles Webber. Therefore, it may contain some differences, which cause that the results may be not fully comparable with the software of Charles Webber.

The RQA software uses a specific normalization of the distance matrix, whereas the CRP toolbox uses either no normalization or a normalization to a standard deviation of one and a mean of zero. However, with an appropriate data preparation and setting you can get the same results.

For compatibility use a Theiler window of size one and ensure that the data are normalized before by the same distance which is used in the RQA software; e.g. normalize with the maximal phase space diameter. In the CRP toolbox, with the programme pss you can estimate the maximal distance of the phase space, which is reconstructed from the data. This max. distance can be used for the normalization of the data:

 d_max = pss(x,3,5,'euclidean');
x_norm = 100*x/d_max;


Applying crqa on this normalized data, using a Theiler window of one and don't let crqa normalize the data (which is the defeault setting!), we will get the same RQA measures as the RQA software of Charles Webber does:

RQA = crqa(x_norm, dim, lag, e, [], [], l_min, v_min, 1, ...
'euclidean', 'nonormalize', 'silent')


top

How can I determine the number of diagonal lines in the RP?

In order to get the number of the diagonal lines, first compute a recurrence plot and assign it to some variable, e.g.

X = crp(x,1,1,.2,'non','silence');


(which means, it uses no embedding, a threshold of 0.2, the data will not be normalized before (what is default) and finally it suppresses output of the calculation).

Then simply use the command dl

[m L] = dl(X);


Now you find in m the mean diagonal line length and in L a vector containing the lengths of each found diagonal line. From L you can get the number of lines with lengths larger than e.g. 2 by

N = length(nonzeros(b>2))

top

How can I determine the histogramme of the lengths of diagonal/vertical lines in the RP?

From the recurrence matrix in X the histogramme of the lengths of diagonal/vertical lines can be determined using the commands dl and tt.

For the diagonal lines simply call

[m L] = dl(X);


and for the vertical lines

[m L] = tt(X);


In m is the mean diagonal/vertical line length and in L a vector containing the lengths of each found diagonal/vertical line. From L the histogramme of the line lengths is then

hist(L,[1:max(L)])

or
H = hist(L,[1:max(L)]);


Please note that this histogramme contains the number of lines which are exactly of length $$l$$.

top

How can I calculate the entropy of the data by means of RQA?

Using RQA it is not possible to calculate the entropy of the data. The entropy measure provides in RQA measures the Shannon entropy of the distribution of the lengths of the diagonal lines in an recurrence plot.

However, using the CRP toolbox you can apply the command entropy to a histogram of the data. This provides you a simple entropy estimation of the data.

top

How can I calculate the RQA trend measure?

There are different ways to calculate the RQA measure TREND using the CRP toolbox. I have not included this measure in the toolbox for several reasons. The main reason is that, in my opinion, this measure depends much stronger than the other measures on the chosen settings. This measure is a bit critical as it can give you results which may be misleading or even contradicting. Therefore, I do not support the distribution and application of this measure by including it in my toolbox.

If you really like to calculate, you can use either CRQAD and (after some cutting of the result vector) you calculate the slope. Or, even more simple, just use the RP-matrix and apply some MATLAB functions. In the following you can find an example where TREND is calculated in sliding windows:

% x = time series
N = 1000;
x = rand(N,1);

m = 3; % embedding dimension
t = 1; % embedding delay
e = .3; % recurrence threshold
w = 200; % window size
ws = 50; % window step

timescale = 1:ws:N - w + 1;

trend = zeros(length(timescale),1);
for i = timescale
X = crp(x(i:i+w-1),m,t,e,'max','non','sil'); % calculate RP
N_X = size(X,2); % size of RP (i.e. number of columns)
T_ = zeros(N_X-1,1);
% count the number of recurrence points in diagonal k
for k = 1:N_X-1
T_(k) = nnz(diag(X,k)) / (N_X-k)*100;
end
Ntau = N_X - 1 - round(0.1*N_X); % last 10% of the RP will be skipped
p = polyfit((2:Ntau+1)',T_(1:Ntau),1); % slope
trend(i) = 1000 * p(1); % Webber's definition includes factor 1000
end

plot(timescale, trend(timescale))


The result is the TREND measure as calculated in Chuck Webber's RQA software.

top

How the RQA measures in the plot of the GUI of the function crqa are aligned to the time scale?

The window of length $$w$$ is applied on the data and not on the RP, i.e., the RP will have smaller size than the window, thus $$w-(m-1)\tau$$. If we consider the data window to time $$i \ldots i+w$$, the corresponding RQA measures are assigned at time $$i$$. Therefore, if you see a beginning of a transition in the plot of the RQA measures at time $$i$$, this transition will probably happen at time $$i+w-(m-1)\tau$$.

top

What is the difference in the output of the command crqad denoted as, e.g., Y.RRp and Y.RRm?

In order to study also anti-correlation, the second time series is additionally multiplied by −1:

$$x(t) \rightarrow x(t)\\ y(t) \rightarrow -y(t)$$

The results of the diagonal-wise computed RQA measures from the CRP are denoted as RRp for the normally constructed CRP and RRm for the CRP based on the "negative" second time series, corresponding to $$RR_+$$ and $$RR_-$$, resp., as introduced in

• N. Marwan, J. Kurths: Nonlinear analysis of bivariate data with cross recurrence plots, Physics Letters A, 302(5–6), 299–307 (2002). DOI:10.1016/S0375-9601(02)01170-2

top

Why the RQA measures contain so many zeros when using a window step size larger than one?

This is just to that we can directly apply the time scale of the original data to the RQA measures. Assume we use a window step size of 10. Then we need to plot only every 10th value:

ws = 10;
Y = crqa(rand(1000,1), 1, 1, 0.1, 50, ws, 'silent');
plot(Y(1:ws:end, 1)


top

How can I predict states using the CRP toolbox like the VRA software of Eugene Kononov?

Using the CRP toolbox it is not possible to predict states. Prediction is not a part of a recurrence analysis. Although Eugene Kononov provides an additional forecast model in his VRA software, which is based on recurrence, we do not intent to include a similar model in the CRP toolbox.

top