H.1 Assignment: Cross Recurrence Quantification Analysis
- Create two variables for CRQA analysis:
y1 <- sin(1:900*2*pi/67) y2 <- sin(.01*(1:900*2*pi/67)^2)
You have just created two sine waves. We’ll examine if and how they are coupled in a shared phase space. As a first step plot them.
Find an embedding delay (using mutual information) and an embedding dimension (if you calculate an embedding dimension using package
fractalfor each signal seperately, as a rule of thumb use the highest embedding dimension you find in further analyses).
# General settings for `crqa()` par0 <- list(rescale = 0, normalize = 0, mindiagline = 2, minvertline = 2, tw = 0, whiteline = FALSE, recpt = FALSE, side = "both", checkl = list(do = FALSE, thrshd = 3, datatype = "categorical",pad = TRUE) ) # Settings for `optimizeParam()` par <- list(lgM = 20, steps = seq(1, 6, 1), radiusspan = 100, radiussample = 40, normalize = par0$normalize, rescale = par0$rescale, mindiagline = par0$mindiagline, minvertline = par0$minvertline, tw = par0$tw, whiteline = par0$whiteline, recpt = par0$recpt, fnnpercent = 10, typeami = "mindip")
- We can now create a cross recurrence matrix. Fill in the values you decided on. You can choose a radius automatically, look in the
Note: There is no rescaling of data, the sines were created in the same range. You can plot a matrix using
image(). You could also check package
nonlinearTseries. If you sourced the
- Get the optimal parameters using a radius which will give us 2%-5% recurrent points.
(ans <- optimizeParam(ts1 = y1, ts2 = y2, par = par, min.rec = 2, max.rec = 5))
$radius  0.2128611 $emddim  2 $delay  15
- Run the CRQA.
crqaOutput <- crqa(ts1 = y1, ts2 = y2, delay = ans$delay, embed = ans$emddim, radius = ans$radius, normalize = par0$normalize, rescale = par0$rescale, mindiagline = par0$mindiagline, minvertline = par0$minvertline, tw = par0$tw, whiteline = par0$whiteline, recpt = par0$recpt, side = par0$side, checkl = par0$checkl )
- Can you understand what is going on? Explain the the lack of recurrent points at the beginning of the time series.
Examine the synchronisation under the diagonal LOS. Look in the manual of
crqaor Coco & Dale (2014).
Make a plot of the diagonal profie. How far is the peak in RR removes from 0 (Line of Synchronisation)?
This is based on the timeseries (without embedding)
win <- 50 (res <- drpdfromts(y1, y2, ws = win, datatype = 'continuous', radius = ans$radius))[2:3]
$maxrec  0.3080357 $maxlag  55
plot(-win:win,res$profile,type = "l", lwd = 5, xlab = "Delays", ylab = "Recurrence") abline(v=0)
To get the diagonal profile from the recurrence plot, use
dprofile <- spdiags(crqaOutput$RP) plot(-win:win,colMeans(dprofile$B[, between(dprofile$d, -win,win)]), type="l", lwd = 5, xlab = "Delays", ylab = "Recurrence") abline(v=0)
- Perform the same steps with a shuffled version (or use surrogate analysis!) of the data of timeseries \(y1\). You can use the embedding parameters you found earlier.