Wavelet Transform & Scattering Transform

--- title: Wavelet Transform & Scattering Transform tags: Practical description: --- # Wavelet Transformation & Scattering Transformation [ToC] ## Wavelet Transform (WT) ### Introduction to Continuous Wavelet Transform (CWT) The basic idea in wavelet analysis is to start with a function $\psi(t)$, known as the mother wavelet. Then, a dilated wavelet is \begin{equation} \psi_{\lambda}(t)=\frac{1}{\lambda}\psi\left(\frac{t}{\lambda}\right) \end{equation} Note that if $\lambda$ is larger then the oscillation frequency is lower. Consider the following function \begin{equation*} \psi_{a, b}(t)=\frac{1}{a} \psi\left(\frac{t-b}{a}\right) \end{equation*} The parameters $a$ and $b$ are responsible for **scaling** and **translating** the function. Hence, we could use above wavelet to plot **scalogram**, with scaling axis and translating axis by \begin{equation} X(a, b)=\int^\infty_{-\infty}f(t)\psi_{a, b}(t)^\ast dt \end{equation} :::warning The higher scale factor is, the lower frequency will be captured. Hence, CWT($a$, $b$) could be viewed as capturing the (frequency, time) at $(\frac{\omega_0}{a},b)$, where $\omega_0$ is origninal frequency of mother wavelet and $\psi$ center at $t=0$. ::: :::danger Sometimes, if define $\psi_{a, b}(t)=\frac{1}{\sqrt{a}} \psi\left(\frac{t-b}{a}\right)$, then there is an factor $\frac{1}{a}$ in the above transformation. ::: * Refer to [S. L. Brunton & J. Nathan Kutz ](http://databookuw.com/databook.pdf) and their [video ](https://www.youtube.com/c/Eigensteve/playlists). * Refer to [S. Qian & D. Chen ](https://www.amazon.com/Joint-Time-Frequency-Analysis-Method-Application/dp/0132543842). ### Focus on Mother Wavelet A complex Morlet wavelet $\psi$ can be defined as the product of a complex sine wave and a Gaussian window \begin{equation} \psi(t)=e^{2i\pi f_0t}e^{-\frac{t^2}{2\sigma_t}}=e^{i\omega_0}e^{-\frac{t^2}{2\sigma_t}} \end{equation} where $f_0$ is frquency and $\sigma_t$ is the bandwidth of the Gaussian in time domain, Actually we have 3 way to determine the Morlet wavelet 1. Gaussian bandwidth ($\sigma_t$) in time domain 2. Gaussian bandwidth ($\sigma_\omega$) in frequency domain 3. number of cycle in bandwidth in time domain :::warning In fact, $\sigma_t$ and $\sigma_\omega$ is satisfied the uncertainty principle. ::: Moreover, the number of cycle ($n$) and $\sigma_t$ follow the following relation \begin{equation} \sigma_t = \frac{n}{2\pi f} \end{equation} On the other hand, Refer to [M. X Cohen (2018) ](https://www.researchgate.net/publication/327144906_A_better_way_to_define_and_describe_Morlet_wavelets_for_time-frequency_analysis) or refer to their [video ](https://www.youtube.com/watch?v=4TTpwIZrUAo&list=PLn0OLiymPak2G__qvavn3T8k7R8ssKxVr&index=3&ab_channel=MikeXCohen). For instance, let $f_0=10$, determine the wavelet by bandwidth in frequency domain. ![](https://i.imgur.com/lE2CIYG.png) :::warning The time frequency trade-off can be seen. ::: ### Implement by bandpass filter Now, apply Fourier transformation on wavelet, say $\hat{\psi}(\omega)$, where $\omega=2\pi f$. Clearly, $\hat{\psi}(\omega)$ is a Gaussian pulse with bandwidth $\sigma_\omega$ and center at $\omega$. It could be viewed as bandpass filter. Now, implement CWT as bandpass filter. \begin{equation} X(a,b)=\frac{1}{2\pi}\int^\infty_{-\infty}\hat{f}(\omega)\hat{\psi}(\omega)e^{i\omega b}dt \end{equation} Refer to [MATLAB help center ](https://www.mathworks.com/help/wavelet/gs/continuous-wavelet-transform-as-a-bandpass-filter.html). ```graphviz digraph graphname{ rankdir="LR"; T [label="signal in TD"] P1 [label="singal in TF"] P2 [label="signal in FD"] Q21 [label="signal in FD"] T->P1 [label="Conv", fontcolor=red] T->P2 [label="FFT", fontcolor=darkgreen] P2->Q21 [label="Mutiply", fontcolor=red] Q21->P1 [label="IFFT", fontcolor=darkgreen] } ``` Please refer to [M. X Cohen video ](https://youtu.be/4TTpwIZrUAo?t=502). ### DWT & Resolution Since the uncertainty principle, the resolution should be as following figure. :::warning First, focus at time domain. If low frequency, the bandwidth $\sigma_\omega$ in frequency domain is small. Hence, the bandwidth $\sigma_t$ in time domain is large, so the time resolution at low frequency is small. And vice versa. ::: ![](https://i.imgur.com/7d3COPP.png) :::warning Second, focus at frequency domain. If high frequency, the bandwidth $\sigma_\omega$ in frequency domain is large. Hence the frequency resolution at high frequency is small. And vice versa. ::: ![](https://i.imgur.com/UhuZAm3.png) According to above figure, we could see that ```graphviz digraph graphname{ rankdir="LR"; T [label="signal"] P1 [label="HF (8Hz)"] P2 [label="LF"] Q21 [label="HF (4Hz)"] Q22 [label="LF"] R221 [label="HF (2Hz)"] R222 [label="LF"] T->P1 [label="WT", fontcolor=red] T->P2 P2->Q21 [label="WT", fontcolor=red] P2->Q22 Q22->R221 [label="WT", fontcolor=red] Q22->R222 } ``` On the other hand, the time frequency trade-off can be seen as following figure. ![](https://i.imgur.com/1xnMXQN.png) According to above figure, we could let $a=2^j$ and $b=k2^j$. The discrete wavelet transformation is as \begin{equation} \psi_{j, k}(t)=\frac{1}{2^j} \psi\left(\frac{t-k 2^j}{2^j}\right) \end{equation} :::danger Note that if $a$ is large then frequency is small. Meanwhile, the space between two time grids is large. ::: ### Octacve In the above equation $j$ might not be integer. We could divide ''octave'' to be finer. The parameter $Q$ is often referred to as the number of ''voices per octave''. \begin{equation} \psi_{j, k}(t)=\frac{1}{2^{j/Q}} \psi\left(\frac{t-k 2^{j/Q}}{2^{j/Q}}\right) \end{equation} where let **new** $\psi(t)$ satisfy $\hat{\psi}(\omega)$ center at frequency $\omega_0$ and with bandwidth $\frac{\sigma_\omega}{Q}$. For instance, if $Q=1$, and $\sigma_\omega=\frac{\omega_0}{2}$, ![](https://i.imgur.com/8Vw3PZ3.png) If $Q=2$, the bandwidth is half as long as original bandwidth. Moreover, it requires two step to reach next octave. ![](https://i.imgur.com/gA5z2Ta.png) The reason $Q$ is referred to as the number of voices per octave is because increasing the scale by an octave (a doubling) requires $Q$ intermediate scales. Refer to [MATLAB help center ](https://www.mathworks.com/help/wavelet/gs/continuous-and-discrete-wavelet-transforms.html). ## Implement Wavelet Transform ### Wavelet Take $Q=2$. Let mother wavelet with $f_0=1$ (Hz) and bandwidth $\sigma_\omega=1/4=\frac{f_0}{2Q}$ in frequency domain. Let wavelet center at $\Lambda=\{2^{j/Q}\}_{j=-4}^{12}$. The following shows that $j$ is even case. ![](https://i.imgur.com/YAjb3Q6.png) ### Convolution with wavelet The first is the original signal.Each row shows that signal convolution with different wavelet.  ![](https://i.imgur.com/OA3GmFR.png) :::warning Note that resolution is different. ::: ### Scalogram Use interpolation to plot the scalogram. ![](https://i.imgur.com/Ukksne5.png)  ## Try MATLAB Built-in Functions Use `cwt`. ![](https://i.imgur.com/oqzlwlK.png) ## Comparison with Synchrosqueezing Transforms (SST) ![](https://i.imgur.com/bBtvZbc.png) ## Scattering Transformation (ST) ### Main idea Extend the graph mentioned. The blue nodes are created by the first WT and the purple nodes are created by the second WT. ```graphviz digraph graphname{ rankdir="LR"; T [label="signal"] P1 [label="HF (8Hz)", color=blue] P2 [label="LF"] Q21 [label="HF (4Hz)", color=blue] Q22 [label="LF"] R221 [label="HF (2Hz)", color=blue] R222 [label="LF"] Q11 [label="HF (4Hz)", color=purple] Q12 [label="LF"] R211 [label="HF (2Hz)", color=purple] R212 [label="LF"] R121 [label="HF (2Hz)", color=purple] R122 [label="LF"] R111 [label="HF (2Hz)"] R112 [label="LF"] T->P1 [label="WT", fontcolor=red] T->P2 P2->Q21 [label="WT", fontcolor=red] P2->Q22 Q22->R221 [label="WT", fontcolor=red] Q22->R222 P1->Q11 [label="WT-2nd-1", fontcolor=blue] P1->Q12 Q21->R211 [label="WT-2nd-2", fontcolor=blue] Q21->R212 Q12->R121 [label="WT-2nd-1", fontcolor=blue] Q12->R122 Q11->R111 [label="WS-3nd", fontcolor=purple] Q11->R112 } ``` :::warning Focus on second layer (blue). There are 3 nodes, which are created by the first wavelet transformation. **However, each node contains different information.** Hence, they convole with different number of wavelet. ::: ## Reference 1. MATLAB Help Center, *Continuous and Discrete Wavelet Transforms*. [](https://www.mathworks.com/help/wavelet/gs/continuous-and-discrete-wavelet-transforms.html) 2. MATLAB Help Center, *Continuous Wavelet Transform as a Bandpass Filter*. [](https://www.mathworks.com/help/wavelet/gs/continuous-wavelet-transform-as-a-bandpass-filter.html) 3. Chun-Lin Liu, *A Tutorial of the Wavelet Transform*, NTU, (2010). [](http://disp.ee.ntu.edu.tw/tutorial/WaveletTutorial.pdf) 4. Steven L. Brunton & J. Nathan Kutz, *Data Driven Science & Engineering: Machine Learning, Dynamical Systems, and Control*, UW, (2017). [](http://databookuw.com/databook.pdf) 5. Steven L. Brunton, Youtube Channels. [](https://www.youtube.com/channel/UCm5mt-A4w61lknZ9lCsZtBw) 6. Kevin E Alexander, *Time Frequency Analysis Plugin for EEGLAB*, Air Force Research Laboratory. [](https://osf.io/sxhpy/) 7. Kevin E Alexander, Youtube Channels. [](https://www.youtube.com/user/0792kevin) 8. Jian-Jiun Ding, *Lecture Note of Time-Frequency Analysis and Wavelet Transform*, NTU, (2021). [](http://djj.ee.ntu.edu.tw/TFW.htm) 9. Chia-Feng Lu, *Lecture Note of MATLAB Programming for Medical Signal Analysis*, NYCU, (2021). [](https://www.ym.edu.tw/~cflu/CFLu_course_matlabsig.html) 10. Rami Khushaba, Youtube Channels, USYD. [](https://www.youtube.com/user/BarwariNSW/featured) 11. Shie Qian, Dapang Chen, *Joint Time-Frequency Analysis: Method and Application*. [](https://www.amazon.com/Joint-Time-Frequency-Analysis-Method-Application/dp/0132543842) 12. Mike X Cohen, Youtube Channel, Radboud University Medical Center. [](https://www.youtube.com/channel/UCUR_LsXk7IYyueSnXcNextQ/playlists) 13. Ingrid Daubechies, *Ten Lectures on Wavelets*. [](https://jqichina.files.wordpress.com/2012/02/ten-lectures-of-waveletsefbc88e5b08fe6b3a2e58d81e8aeb2efbc891.pdf)