10  Wavelet Transform

With Wavelets we extend the concept of the Fourier analysis to general orthogonal bases. This extensions is done in such a way that we can do a multi-resolution decomposition and thus partially overcome the uncertainty principal discussed before.

Wavelets are both, local and orthogonal. The whole family of a wavelet are generated by scaling and translating a mother wavelet \(\psi(t)\) as \[ \psi_{a,b}(t) = \frac{1}{\sqrt{a}} \psi\left(\frac{t-b}{a}\right), \] where the parameters \(a\) and \(b\) are responsible for scaling and translating, respectively.

The simplest example is the so called Haar wavelet.

Example 10.1 (The Haar wavelet) The mother wavelet is defined as the step function \[ \psi(t) = \begin{cases} \begin{array}{rl} 1, & \text{for}\, 0 \leq t \le \tfrac12\\ -1, & \text{for}\, \tfrac12 \leq t \le 1\\ 0, & \text{else} \end{array} \end{cases} \]

Show the code for the figure
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
%config InlineBackend.figure_formats = ["svg"]

N = 1000
a, b = -0.25, 1.25
t = np.linspace(a, b, N, endpoint=False)
mother = lambda t: np.heaviside(t, 1) - np.heaviside(t - 1/2, 1) - \
                   (np.heaviside(t - 1/2, 1) - np.heaviside(t - 1, 1))
psi = lambda t, a, b: 1 / np.sqrt(a) * mother((t - b) / a)

plt.figure(0)
plt.plot(t, psi(t, 1, 0))
plt.xlabel("t")
plt.ylabel(r"$\psi_{1,0}$")
plt.gca().set_aspect(0.125)

plt.figure(1)
plt.plot(t, psi(t, 1/2, 0))
plt.xlabel("x")
plt.ylabel(r"$\psi_{\frac{1}{2}, 0}$")
plt.gca().set_aspect(0.125)

plt.figure(2)
plt.plot(t, psi(t, 1/2, 1/2))
plt.xlabel("t")
plt.ylabel(r"$\psi_{\frac{1}{2}, \frac{1}{2}}$")
plt.gca().set_aspect(0.125)

plt.show()
(a) Scaling 1 and translation 0.
(b) Scaling 1/2 and translation 0.
(c) Scaling 1/2 and translation 1/2.
Figure 10.1: Haar wavelets for the first two levels of multi resolution.

Note that the Haar wavelets are orthogonal and provide a hierarchical basis for a signal.

Example 10.2 (The Maxican hat wavelet) The mother wavelet is defined as the negative normalized second derivative of a Gaussian function, \[ \psi(t) = (1 - t)^2\, \mathrm{e}^{-\tfrac{t}{2}} \]

Show the code for the figure
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
%config InlineBackend.figure_formats = ["svg"]

N = 1000
a, b = -5, 5
t = np.linspace(a, b, N, endpoint=False)
mother = lambda t, d: (1 - np.pow(t, 2)) * np.exp(-1/2 * np.pow(t, 2))
psi = lambda t, a, b: 1 / np.sqrt(a) * mother((t - b) / a, 2)

plt.figure(0)
plt.plot(t, psi(t, 1, 0))
plt.xlabel("t")
plt.ylabel(r"$\psi_{1,0}$")
plt.gca().set_aspect(1.5)

plt.figure(1)
plt.plot(t, psi(t, 1/2, 0))
plt.xlabel("x")
plt.ylabel(r"$\psi_{\frac{1}{2}, 0}$")
plt.gca().set_aspect(1.5)

plt.figure(2)
plt.plot(t, psi(t, 1/2, 1/2))
plt.xlabel("t")
plt.ylabel(r"$\psi_{\frac{1}{2}, \frac{1}{2}}$")
plt.gca().set_aspect(1.5)

plt.show()
(a) Scaling 1 and translation 0.
(b) Scaling 1/2 and translation 0.
(c) Scaling 1/2 and translation 1/2.
Figure 10.2: Mexican hat wavelets for the first two levels of multi resolution.

If we have a wavelet \(\psi\), we can generate a new wavelet by convolution \[ \psi \ast \phi, \] for a bounded and integrable function \(\phi\).

Definition 10.1 (Continuous Wavelet Transform) The continuous wavelet transform is given by \[ \mathcal{W}_\psi\{f\}(a, b) = \langle f, \psi_{a, b} \rangle = \int_{-\infty}^\infty f(t)\overline{\psi}_{a,b}\, \mathrm{d}t, \] this is only true for a bounded wavelet \(\psi\) \[ C_\psi = \int_{-\infty}^\infty \frac{|\hat{\psi}(\tau)|^2}{|\tau|}\, \mathrm{d}\tau, \] i.e. a wavelet with \(C_\psi < \infty\). In this case also the inverse transform exists and is defined as \[ f(t) = \frac{1}{C_\psi} \int_{-\infty}^\infty \int_{-\infty}^\infty\mathcal{W}_\psi\{f\}(a, b)\psi_{a, b}(t)\frac{1}{a^2} \, \mathrm{d}a\, \mathrm{d}b, \]

From the continuous wavelet transform we go on to the discrete wavelet transform as similar for the Fourier transforms we have seen we will hardly ever have the entire signal at hand for the transformation.

Definition 10.2 (Discrete Wavelet Transform) The discrete wavelet transform is given by \[ \mathcal{W}_\psi\{f\}(j, k) = \langle f, \psi_{j, k} \rangle = \int_{-\infty}^\infty f(t)\overline{\psi}_{j, k}\, \mathrm{d}t, \] where \(\psi_{j, k}\) is a discrete family of wavelets \[ \psi_{j, k}(t) = \frac{1}{a^j}\psi\left(\frac{t - k b}{a^j}\right). \] The inverse is than given by \[ f(t) = \sum_{j,k = -\infty}^\infty \mathcal{W}_\psi\{f\}(j, k) \psi_{j, k}(t). \] Which is nothing else than expressing the function in the wavelet family. If this family of wavelets is orthogonal (as e.g. the Haar wavelet) it is possible to expand the function \(f\) uniquely as they form a basis.

Note

There also exists a fast wavelet transform that reduces the computational complexity from \(\mathcal{O}(N\log N)\) to \(\mathcal{O}(N)\) by cleverly reusing parts of the inner product computation.

Example 10.3 (Signal analysis with the Haar wavelet) To start and get an idea how the analysis works we use an instructive example. We have a piecewise constant function as \(v = [3, 1, 0, 4, 0, 6, 9, 9]^{\mathsf{T}}\) \[ f(x) = \sum_{i=1}^{8} v_i \chi_{[i-1, i)}, \] where \(\chi_{[c, d)}\) is the indicator function that is 1 on the interval \([c, d)\) and zero everywhere else.

Now let us proceed through the levels of the transform to see what is happening, we select \(a=2\) and \(b=1\):

\[ \mathcal{W}\{f\}(1, :) = \left( \begin{array}{cccccccc} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0&0&0&0&0\\ 0&0&\frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0&0&0\\ 0&0&0&0&\frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0\\ 0&0&0&0&0&0&\frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}&0&0&0&0&0&0\\ 0&0&\frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}&0&0&0&0\\ 0&0&0&0&\frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}&0&0\\ 0&0&0&0&0&0&\frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}\\ \end{array} \right) \left[ \begin{array}{c} 3\\ 1\\ 0\\ 4\\ 0\\ 6\\ 9\\ 9 \end{array} \right]=\left[ \begin{array}{c} 2\sqrt{2}\\ 2\sqrt{2}\\ 3\sqrt{2}\\ 9\sqrt{2}\\ \sqrt{2}\\ -2\sqrt{2}\\ -3\sqrt{2}\\ 0 \end{array} \right] \]

We then apply the next level just for the first half \[ \mathcal{W}\{f\}(2, :) = \left( \begin{array}{cccccccc} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0&0&0&0&0\\ 0&0&\frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0&0&0\\ \frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}&0&0&0&0&0&0\\ 0&0&\frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}&0&0&0&0\\ 0&0&0&0&1&0&0&0\\ 0&0&0&0&0&1&0&0\\ 0&0&0&0&0&0&1&0\\ 0&0&0&0&0&0&0&1 \end{array} \right) \left[ \begin{array}{c} 2\sqrt{2}\\ 2\sqrt{2}\\ 3\sqrt{2}\\ 9\sqrt{2}\\ \sqrt{2}\\ -2\sqrt{2}\\ -3\sqrt{2}\\ 0 \end{array} \right]=\left[ \begin{array}{c} 4\\ 12\\ 0\\ -6\\ \sqrt{2}\\ -2\sqrt{2}\\ -3\sqrt{2}\\ 0 \end{array} \right] \] and our final step \[ \mathcal{W}\{f\}(3, :) = \left( \begin{array}{cccccccc} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0&0&0&0&0\\ \frac{1}{\sqrt{2}}&\frac{-1}{\sqrt{2}}&0&0&0&0&0&0\\ 0&0&1&0&0&0&0&0\\ 0&0&0&1&0&0&0&0\\ 0&0&0&0&1&0&0&0\\ 0&0&0&0&0&1&0&0\\ 0&0&0&0&0&0&1&0\\ 0&0&0&0&0&0&0&1 \end{array} \right) \left[ \begin{array}{c} 4\\ 12\\ 0\\ -6\\ \sqrt{2}\\ -2\sqrt{2}\\ -3\sqrt{2}\\ 0 \end{array} \right]=\left[ \begin{array}{c} 8\sqrt{2}\\ -4\sqrt{2}\\ 0\\ -6\\ \sqrt{2}\\ -2\sqrt{2}\\ -3\sqrt{2}\\ 0 \end{array} \right] \]

Show the code for the figure
import numpy as np
import matplotlib.pyplot as plt
import pywt
%config InlineBackend.figure_formats = ["svg"]


N = 1000
a, b = 0, 8
t = np.linspace(a, b, N, endpoint=False)
v = np.array([3, 1, 0, 4, 0, 6, 9, 9])
chi = lambda t, a, b: np.heaviside(t - a, 1) - np.heaviside(t - b, 1)
mother = lambda t: chi(t, 0, 1/2) - chi(1/2, 1)
psi = lambda t, a, b: 1 / np.sqrt(a) * mother((t - b) / a, 2)

def f(t, v, spread=1):
    y = np.zeros(t.shape)
    for i, x in enumerate(v):
        y += x * chi(t, i * spread, (i + 1) * spread)
    return y

X = np.zeros((4, len(v)))
for i in range(0, 4):
    x = pywt.wavedec(v, wavelet="Haar", level=i)
    X[i, :] = np.concatenate(x) 

plt.figure(0)
plt.plot(t, f(t, v), label="Signal")
#plt.plot(t, f(t, X[0, :], 1))
plt.plot(t, f(t, X[1, 0:5], 2), label="Level 1")
plt.plot(t, f(t, X[2, 0:3], 4), label="Level 2")
plt.plot(t, f(t, X[3, 0:1], 8), label="Level 3")
plt.xlabel("t")
plt.legend()
plt.show()
Figure 10.3: Haar wavelets for the first three levels of multi resolution analysis.