High dimensional chaos (HDS).
Welcome to the Web Site dedicated to the algorithm
for estimation of the correlation dimension d of High Dimensional Signals.
The new version has been uploaded on 29 Dec 2018. Some new details have been added and the algorithm has been accelerated.
What has been modified in the version-2 of the HDS toolkit? New algorithm details are marked with (new) ....
The modifications marked with (new) were developed after publishing the papers below and are not described in the cited papers.
1. K. Michalak, (2011) Modifications of the Takens-Ellner algorithm for medium- and high-dimensional signals, Physical Review E, Vol. 83(2):026206,
2. K.P. Michalak. How to estimate the correlation dimension of high-dimensional signals? AIP Chaos, 24(2):033118, 2014.
3. K.P. Michalak. Estimating correlation dimension of high-dimensional signals - quick algorithm. AIP Advances 8, 105201 (2018); doi: 10.1063/1.5013255
The notation being used in the functions dhds.m and dhdsall.m corresponds to the notation being used in the above mentioned articles.
DOWNLOAD HDS-TOOLKIT-v2 (50MB) (updated 29.12.2018)
This is the beta version of the toolkit.
If something does not go, please write me the email with the description of your problem. I try help.
Download the former version 1 if necessary...
There are 2 functions dedicated for estimating correlation dimension d
dhds.m - the short version of algorithm is implemented
dhdsall.m - the full version of algorithm is implemented
The general call of dhds.m function is:
dhds (series, quick precise, MaxIterationTime, Da, DeflErr, norm, NW Wmax WV, Pmin, dcW, fitting, figures)
However, you can use the simplest call: d = dhds ( series )
Use empty string: '' if you want use the default values of subsequent input parameters.
Let us analyze the input parameters:
(1) series - the signal being analyzed,
use empty string or nothing for looking the exemplary calculation of the 3xLorenz signal having d = about 6 (new).
(2) quick precise_surrog: this parameter joins 2 parapmeters: speed and surrogate number .(new)
In general, please use
the number in the format ssp
: ss -
surrogates number, p = 1, 2
Use 1 for preliminary estimation and 3 for final estimation, if necessary. Change other parameters if the results are not satisfactory.
- 1 or
'quick' - quick, low precision, no surrogates
- e.g 11
- quick, 1 surrogate
It is possible, also, to change the
HEADER's variable 'UserAction'. If 'UserAction'=0 than the calculations
are repeated using increased Deflection Error e. If 'UserAction'=1, the algorithm
stops and asks for the user's decision: 'Follow' using current e'
'Break and repeat using increased e'.
(4) Da: assumed nbr of distances being shorter than the estimated distance dmax. dmax is estimated based on assumed deflection error e and used embedding dimension m. When omitted, default value of Da depends on the used quick_precise_surrog parameter: Da = 1000 or 3000 or 10000.
It must be pointed, however, that to small e may shift the results into `noise` region of d=fn(Pmax) relation.
You can use higher initial values e.g. 5-10%
- if the signal is significantly charged with noise,
- if the expected complexity d is very high
- if the signal is short
- if you want obtain robust results in very quick time.
You may use lower values if the complexity of the signal is rather low, the signal is noise-free (generated by some differential equations) and you would like reduce the error connected with e-correction. High values accelerate strongly the calculations for small Ws.
The calculation time for large Ws (large d estimated) should be rather regulated by using MaxIterationTime
The 2nd element of this input parameter DeflErr(2) = ecut represents the deflection error for which the dcut is estimated. It must be equal or greater than DeflErr(1)=e. dcut determines the cut-off distance and, thus, the length of the actual dmax vector in the memory space and the duration of its sorting, after Da is reached. The distances shorter than dcut and longer than dmax are the reserve if increase in e is necessary for the given W-iteration. The iteration has not to be repeated. It influences also the Figures being printed as a output of the function. If you use higher value, the relation d = fn (Pmax) is broader. The acceptable values lies between 1 and 30. The unit is [%]. By default, it is 5%.
Both these parameters are
automatically increased by 1 if the assumed number of distances
Da being shorter
than dmax is not reached in the assumed
(7) NW_Wmax_WV: Here, you can set the vector WV with your Window Width values used for calculations. You can just set the vector WV with the assumed W values. E.g. you want to use the same W values for different similar signals.
If the number of elements is less than 4
the algorithm finds default values depending on the number of elements:
Wmax - maximal value of the
vector WV, by default 10t (autocorelation
time of the signal).
- if the parameter is 2-dimensional matrix (2 x NW) - the second row is treated as default embedding dimension m for W's defined in the first row (new)
- if the parameter is 3-dimensional matrix (3 x NW)- the second row - as above - embedding dimension, third row is treated as J for W's defined in the first row and m defined in the second row. (new)
- if [0,Wmax] - dw's are calculated for increasing W's
calculated using a rule:
sets 40 values W from 0.1 to 10 due to estimated autocorrelation time of
the signal, unequally distributed with the maximum density for about 2.5
The interpretation of this parameter depends on its value:
- If 0 < Pmin < 0.1 - it is
treated as Pmin. You can set this value if you have observed the begin
of the influence of noise on the graph d=fn(Pmax,W) in the preliminary
run of the algorithm.
Look for the variable PminCut=1.5 in the HEADER of the function to modify this range.
E.g. PminCut=1.2 denotes the range dinit/1.2 < dinit < dinit* 1.2
It is observed that the relation d=fn(Pmax)
possesses lacunarities in the plateau region. Thus, in order to increase
the precision, a number of points from this relation should be taken to
estimate dw as a mean over these points. The parameter
the the choice of the range of indices Imaxdown
- Imaxup which are used to estimate this mean.
First: Imaxopt is estimated based on dinit, e and m. Iminopt is estimated based on Pmin.
(additional small modifications may be
added to these rules for small Imaxopt if necessary)
(10) fitting_method: This parameter determines the method of polynomial
fitting applied to the relation d=fn(W).
Here, you can define other functions
polynomials to find optimal fitting. The procedure works on symbolic
expressions and solving differential equations, thus be careful with defining
variables to be optimized (especially, do not use the variables in the
denominanitve - symbolic solver fails in this case).
(10) figures :
(12) signal_number : important only for figure names to to create unique figure names when the function is repeated in the loop. This number is added to figure names to distinguish figure names while calling to function in the loop for subsequent signals. Other case - figure are over-saved for subsequent calls of the function.
(13) fighandleset: vector with fighandles to use while printing figures.
when ommitted or 0 or -1 or NaN - new figures are created for every function loop. However, when using signal_number - they are saved with unique names.
(1:3) handles for 3 main figures printed and refreshed after every W-iteration
(4) handle for figure showing the difference between d for original and shuffled signals , only if the number of surrogates defined in parameter(2): quick precise_surrog is higher than 0
(5-10) - dedicated for future use
(11-..) handles for figures with polynomial fittig , odd number - without e-correction, even number - with e-correction
(11-12) first defined in (10)fitting_method polynomial degree or for the optimal_fitting procedure fitting_method =0
(13-14) second defined in (10)fitting_method polynomial degree
(15-16) third defined in (10)fitting_method polynomial degree
(14) -EachToEach (new) - defines the algorithm for the choice of pairs of points in m-dimensional space to calculate distances...
- if 1 - the approach `Each-To-Each` from the A matrix is used to calculate distances. However, the algorithm stops after reaching Da distances shorter than dmax. In order to avoid biases, the order of basal points for comparing with other ones are chosen randomly.
Based on the former W-iteration, the automatically estimated parameter Jreduction lowers the preliminary estimated J if A matrix is expected to be too small to get sufficient number of distances Da being shorter than dmax.
(Lower J -> larger A matrix -> more distances to calculate using EachToEach approach)
- if 0 - (default) repeating the Random Joining Procedure is used, as described in [1,2,3]. This approach may lead to repeating the pairs in subsequent RJP repetitions. However, if the number of RJP repetitions becomes too large, the algorithm jumps automatically to EachToEach approach in next W-iterations. If EachToEach becomes insuficient - J is lowered to get bigger A matrix with m-dimensional points.
(new): Due to the repeatition of the algorithm for subsequent surrogates, the necessity occured to store the results for all surrogate repetitions. Thus, the index s is added in the output paramters.
s = 1 - the reuslts for original signal,
s = 2 .. NSurr+1 - the results for 1..NSurr surrogates, NSurr - the number of surrogates defined in quick_precise_surrog variable.
Used for printing d=fn(Pmax,W)
Analyzed signal - is the sum of 3 Lorenz signals possessing similar time scales.
The short relation d=fn(Pmax,W) after 40 W-iterations: (Figure 2 in the alorithm)
The exemplary relation d=fn(Pmax) for the last W-iteration (Figure 1 in the algorithm)
The upper relation (with e-correction) reflects the lower one (without e-correction) being increased by about 10% due to the used deflection error e=10% corresponding to such very small Pmax's.
Exemplary 6-degree polynomial fitting applied to dw = fn(W) relation after finishing the calculations for the original signal
Comparison with the surrogates:
The difference between <dsurr> and dorig . The increase in dsurr after shuffling the phases in fft (surrogate signals) points to the chaotic structure which is visible especially in the range of W = 400-800 px.
If you would like to be informed about the new events on this web site, please write the email to the author : email@example.com