Can you perform a wavelet coherence calculation with synchrosqueezing of the signals?

Hi,

I'm analyzig some dynamic time domain data and comparing multiple signals. So far, I have good results using the wavelet transform (cwt), but I also tried the synchrosqueezed transform (wsst) and I'm exhilerated by it's effects on each signal! However, what I really need is the wavelet coherence (wcoherence), and I was wondering if there's a way to integrate the synchrosqueezing into wcoherence so that the two signals wsst'ed before being wcoherence'd.

Answer 1

[Gayatri Rathod]

Hi Vitek,

Yes, it is possible to perform a wavelet coherence calculation with synchrosqueezing of the signals in MATLAB. You can use the Wavelet Toolbox, which provides functions for computing the wavelet transform, wavelet coherence, and synchrosqueezing transform.

Here's an example of how you can compute the wavelet coherence with synchrosqueezing:

% Load the data for the two signals you want to analyze. 
load signal1.mat 
load signal2.mat 

% Compute the synchrosqueezed wavelet transform of each signal using the wsst function.
scales = 1:128; 
freqs = logspace(log10(0.01), log10(1), 128); 
wst1 = wsst(signal1, scales, freqs); 
wst2 = wsst(signal2, scales, freqs); 
 
% Compute the wavelet coherence of the two signals using the wcoherence function, using the synchrosqueezed wavelet transform as input. 
wcoh = wcoherence(wst1, wst2); 

• The resulting wcoh matrix will contain the wavelet coherence values for each frequency and scale.
• It's important to remember that synchrosqueezing could result in certain alterations in the frequency domain., which could affect the coherence calculation. You may want to compare the results with and without synchrosqueezing to ensure that it's appropriate for your specific application.

You can read more about the wsst and wcoherence from the following documentations:  wsst Function, wcoherence Function.

Hope it helps!  

Regards,

Gayatri Rathod

[Vitek Stepien]

Are you sure about that? I read the documentation for both functions pretty thoroughly, and I found no indication that you can use the wavelet transform as input to wcoherence. Even looking at the code for wcoherence, it takes the inputs and performs the wavelet transform on them. Wouldn't the method you're suggesting involve taking the wavelet coherence twice?

Edit: I tried using your code with sample 1D real signals, and as I suspected I got the following error:

Error using wcoherence
Expected X to be real.

Error in wcoherence (line 144)
validateattributes(x,{'single','double'},{'real','finite'},...
    

Inputs to wcoherence have to be 1 D real signals.

The way I was thinking could work is if I made a copy of wcoherence and manualy wrote the synchrosqueezing inside the code, replacing the regular cwt, but I wanted to avoid that route.

Answer 2

[xiaoran Yin]

I'm sure the above code is not feasible because that's what I was thinking at the beginning, now I've modified the wsst by combining it with the settings of the wavelet coherence toolbox so that the results obtained can be generated including AR1 significant tests, COI ranges, etc. (of course these are in the tool above) The results are quite promising, at least up to this step of generating CWT photo.

left is cwt and right is sst-cwt, I purposely did not modify the range of significant regions because it is not appropriate, in my opinion, to perform significance tests on the compressed coefficients.

However, the results are not so promising for the next step, which is the coherence calculation of the two compressed wavelet coefficient matrices.The figure below shows the results of my wavelet coherence of the two time series, at this point it can be seen that compared to the original results, what we want does not seem to be reflected, i.e. the wavelet coefficients get compressed in frequency in the high region, but simply disappear.

My following idea is to perform a similar compression on the wavelet coherence coefficients as I did on the continuous wavelet coefficients, but since I am not a math or signal major, and in fact I am classified as a researcher in the field of geosciences, this process has become increasingly difficult for me, and I would be grateful if anyone could help me with this part.

In addition, the above if there is any reading difficulties that is normal, because this is the translated by software translation, I really get the code to get a headache. ;(

You can contact me through my email:[e-mail] or [e-mail]

@Gayatri Rathod@Vitek Stepien hope it can help you.

[Vitek Stepien]

@xiaoran Yin @Gayatri Rathod I don't know if it helps you, but I figured it out. You basically have to write your own wcoherence function, based on the built-in MATLAB script. I copied all the content of wcoherence to a new function wcoherence_mod, not to mess with the original script, and just dug through the lines of code. You can insert extra parts into the code to replace the cwt's with sst's, so that the final coherence caluclation is based on the synchrosqueezed transforms.

Of course there is plenty of parameters you have to keep tabs on, like the number of voices, number of scales, etc. but I'm attaching my modified function here. To use the SST based version, add input arguments pair 'transformtype','wsst' like this:

[wcoh, wcs, freqs, coi] = wcoherence_mod(in1, in2, Fs, 'transformtype', 'wsst');

P.S. I hope that sharing a modified MATLAB function isn't an infringment of any copyright. If it is, I am unaware of that at the time of posting and didn't intend any harm and will remove the function from this post as soon as I'm notified.

[xiaoran Yin]

@Vitek StepienSo glad this issue is getting attention again.Actually your idea is the same as I answered at the beginning, but after I answered it, I analyzed it in depth later on. WTC's idea can be understood as

1. cwt result of two sequences to get Wx,Wy.

2. xwt analyze it to get Wxy(Wxy=Wx.*conj(Wy);).

3. use the smoothwavelet to get sWxy, sWx, sWy.

4. compute Rsq, Rsq=abs(sWxy). ^2./(sX.*sY);

For the parameters such as octave and so on in this is not essential, and the phase angle is from the angle(Wxy) function which is also the important results. For compressing the Rsq wavelet coherence coefficients of the WTC, what we need to do is not to compress the result of the cwt from wsst. I ran your code, and the results I got are roughly the same from my figure. In my opinion, the core of achieving synchronous compression of WTC is to compress Rsq as cwt is compressed by wsst, i.e., rearrangement of wavelet coherence coefficients as in wsst.

Function wsst in matlab at last provides the function sstalgo to do the rearrangement,which is the algo of wsst (note: some modification is needed in this case if you are using a custom freq). sstalgo has three input variables, cwtcfs,phasetf,gamma.

The first is cwt coefficients, that is, the need to rearrange the object, here can be replaced by Rsq.

The third is the minimum value of cwt coefficients, you can do not care about.

The third is the core, phasetf, the wavelet transform coefficients and their derivatives of the phase difference (phasetf = imag(dcwtcfs./cwtcfs). /(2*pi);).

Undoubtedly Rsq itself doesn't have this imaginary part for us to compute, so what's possible is to compute it from inside the results of the Wxy from XWT, but the results are still unsatisfactory.

I have given up on this part for many reasons. I will upload the code I have modified, but since I am not a professional programmer, the code is very messy. And note that I am not using the wavelet tool that comes with matlab, but the Wavelet coherence toolbox from Grinsted(which is much more friendly for geostudy,GitHub - grinsted/wavelet-coherence: A cross wavelet and wavelet coherence toolbox for MATLAB),

so I hope this helps.

P.S. I hope that sharing a modified MATLAB function isn't an infringment of any copyright. If it is, I am unaware of that at the time of posting and didn't intend any harm and will remove the function from this post as soon as I'm notified( ;)Just like you. )