Subsections

• 7.3.1 A Case Study
• 7.3.2 Continuous Wavelet Transform
  • Resolution of Identity.
  
  
• 7.3.3 Multiresolution Analysis
  • Normalization.
  • Orthogonality.
  
  
• 7.3.4 Haar Wavelets
• 7.3.5 Daubechies' Wavelets

7.3 Wavelets
and Other Multiscale Transforms

Given their recent popularity and clear evidence of wide applicability the most of the space in this chapter is devoted to Wavelet transforms. Statistical multiscale modeling has, in recent decade, become a well established area in both theoretical and applied statistics, with impact to developments in statistical methodology.

Wavelet-based methods are important in statistics in areas such as regression, density and function estimation, factor analysis, modeling and forecasting in time series analysis, in assessing self-similarity and fractality in data, in spatial statistics.

The attention of the statistical community was attracted in late 1980's and early 1990's, when Donoho, Johnstone, and their coauthors demonstrated that wavelet thresholding, a simple denoising procedure, had desirable statistical optimality properties. Since then, wavelets have proved useful in many statistical disciplines, notably in nonparametric statistics and time series analysis. Bayesian concepts and modeling approaches have, more recently, been identified as providing promising contexts for wavelet-based denoising applications.

In addition to replacing traditional orthonormal bases in a variety statistical problems, wavelets brought novel techniques and invigorated some of the existing ones.

7.3.1 A Case Study

We start first with a statistical application of wavelet transforms. This example emphasizes specificity of wavelet-based denoising not shared by standard state-of-art denoising techniques.

A researcher in geology was interested in predicting earthquakes by the level of water in nearby wells. She had a large ( $8192=2^{13}$ measurements) data set of water levels taken every hour in a period of time of about one year in a California well. Here is the description of the problem.

The ability of water wells to act as strain meters has been observed for centuries. The Chinese, for example, have records of water flowing from wells prior to earthquakes. Lab studies indicate that a seismic slip occurs along a fault prior to rupture. Recent work has attempted to quantify this response, in an effort to use water wells as sensitive indicators of volumetric strain. If this is possible, water wells could aid in earthquake prediction by sensing precursory earthquake strain.

We have water level records from six wells in southern California, collected over a six year time span. At least $13$ moderate size earthquakes (magnitude ${4.0}$-${6.0}$) occurred in close proximity to the wells during this time interval. There is a significant amount of noise in the water level record which must first be filtered out. Environmental factors such as earth tides and atmospheric pressure create noise with frequencies ranging from seasonal to semidiurnal. The amount of rainfall also affects the water level, as do surface loading, pumping, recharge (such as an increase in water level due to irrigation), and sonic booms, to name a few. Once the noise is subtracted from the signal, the record can be analyzed for changes in water level, either an increase or a decrease depending upon whether the aquifer is experiencing a tensile or compressional volume strain, just prior to an earthquake.

A plot of the raw data for hourly measurements over one year ( $8192=2^{13}$ observations) is given in Fig. 7.6a, with a close-up in Fig. 7.6b. After applying the wavelet transform and further processing the wavelet coefficients (thresholding), we obtained a fairly clean signal with a big jump at the earthquake time. The wavelet-denoised data are given in Fig. 7.6d. The magnitude of the water level change at the earthquake time did not get distorted in contrast to traditional smoothing techniques. This local adaptivity is a desirable feature of wavelet methods.

For example, Fig. 7.6c, is denoised signal after applying supsmo smoothing procedure. Note that the earthquake jump is smoothed, as well.

Figure 7.6: (a) shows $n=8192$ hourly measurements of the water level for a well in an earthquake zone. Notice the wide range of water levels at the time of an earthquake around $t=417$. (b) focusses on the data around the earthquake time. (c) demonstrates action of a standard smoother supsmo, and (d) gives a wavelet based reconstruction

(a)

(b)

(c)

(d)

7.3.2 Continuous Wavelet Transform

The first theoretical results in wavelets had been concerned with continuous wavelet decompositions of functions and go back to the early 1980s. Papers of Morlet et al. (1982)[19] and Grossmann and Morlet (1984, 1985)[13,14] were among the first on this subject.

Let $\psi_{a,b}(x)$, $a \in {{\mathbb{R}}} \backslash \{0\}, b \in {{\mathbb{R}}}$ be a family of functions defined as translations and re-scales of a single function $\psi(x) \in {{\mathbb{L}}}_2({{\mathbb{R}}})$,

$$\psi_{a,b} (x) = \frac{1}{\sqrt{\vert a\vert}} \psi\left(\frac{x-b}{a}\right){}.$$ (7.16)

Normalization constant $1/\sqrt{\vert a\vert}$ ensures that the norm $\vert\vert\psi_{a,b} (x)\vert\vert$ is independent of $a$ and $b.$ The function $\psi$ (called the wavelet function) is assumed to satisfy the admissibility condition,

$$C_{\psi} = \int_{{{\mathbb{R}}}} \frac{\vert\Psi(\omega)\vert^2}{\vert\omega\vert} \mathrm{d} \omega < \infty{},$$ (7.17)

where $\Psi(\omega)=\int_R \psi(x) \mathrm{e}^{-\mathrm{i} x \omega} \mathrm{d} x$ is the Fourier transform of $\psi(x).$ The admissibility condition (7.17) implies

$$0 = \Psi(0) = \int \psi(x) d\mathrm{d} x{}.$$

Also, if $\int \psi(x) \mathrm{d} x = 0$ and $\int (1+\vert x\vert^{\alpha}) \vert \psi(x)\vert \mathrm{d} x < \infty$ for some $\alpha>0$, then $C_{\psi} < \infty$.

Wavelet functions are usually normalized to ''have unit energy'', i.e., $\vert\vert \psi_{a,b}(x)\vert\vert =1$.

For example, the second derivative of the Gaussian function,

$$\psi(x) = \frac{\mathrm{d}^2}{\mathrm{d} x^2} \left[ - C \mathrm{... ...(1-x^2\right) \mathrm{e}^{-x^2/2}{}, \quad C=\frac{2}{\sqrt{3 \sqrt{\pi}~} }{},$$

is an example of an admissible wavelet, called Mexican Hat or Marr's wavelet, see Fig. 7.7.

Figure 7.7: Mexican hat wavelet (solid line) and its Fourier transform (dashed line)

For any square integrable function $f(x)$, the continuous wavelet transform is defined as a function of two variables

$${\mathcal{CWT}}_f(a,b)= \langle f, \psi_{a,b} \rangle = \int f(x) \overline {\psi_{a,b}(x)} \mathrm{d} x{}.$$

Here the dilation and translation parameters, $a$ and $b$, respectively, vary continuously over ${{\mathbb{R}}} \backslash \{0\} \times {{\mathbb{R}}}$.

Figure 7.8 gives the doppler test function, $f= 1/(t+{0.05}) \sqrt{t (1-t)} \sin(2 \pi \cdot {1.05})$, $0 \leq t \leq 1$, and its continuous wavelet transform. The wavelet used was Mexican Hat. Notice the distribution of ''energy'' in the time/frequency plane in Fig. 7.8b.

Figure 7.8: (a) Doppler signal; (b) Continuous wavelet transform of doppler signal by the Mexican hat wavelet

(a)
(b)

Resolution of Identity.

When the admissibility condition is satisfied, i.e., $C_{\psi} < \infty$, it is possible to find the inverse continuous transform via the relation known as resolution of identity or Calderón's reproducing identity,

$$f(x) = \frac{1}{C_{\psi}} \int_{{{\mathbb{R}}}^2} {\mathcal{CWT}}_f(a,b) \psi_{a,b}(x) \frac{\mathrm{d} a\,\mathrm{d} b}{a^2}{}.$$

The continuous wavelet transform of a function of one variable is a function of two variables. Clearly, the transform is redundant. To ''minimize'' the transform one can select discrete values of $a$ and $b$ and still have a lossless transform. This is achieved by so called critical sampling.

The critical sampling defined by

$$a = 2^{-j}{},\quad b=k 2^{-j}{}, \quad j, k \in {{\mathbb{Z}}}{},$$ (7.18)

will produce the minimal, but complete basis. Any coarser sampling will not produce a unique inverse transform. Moreover under mild conditions on the wavelet function $\psi$, such sampling produces an orthogonal basis $\{\psi_{jk}(x) = 2^{j/2} \psi(2^{j} x -k)$, $j,k \in {{\mathbb{Z}}} \}$. To formally describe properties of minimal and orthogonal wavelet bases a multiresolution formalism is needed.

7.3.3 Multiresolution Analysis

Fundamental for construction of critically sampled orthogonal wavelets is a notion of multiresolution analysis introduced by Mallat (1989a, 1989b)[16,17]. A multiresolution analysis (MRA) is a sequence of closed subspaces $V_n, n \in {{\mathbb{Z}}}$ in ${{\mathbb{L}}}_2({{\mathbb{R}}})$ such that they lie in a containment hierarchy

$$\cdots \subset V_{-2} \subset V_{-1} \subset V_0 \subset V_{1} \subset V_{2} \subset \cdots{}.$$ (7.19)

The nested spaces have an intersection that contains only the zero function and a union that contains all square integrable functions.

$$\cap_n V_j = \{ \boldsymbol{0} \}{}, \quad\overline{\cup_j V_j} = {{\mathbb{L}}}_2({{\mathbb{R}}}){}.$$

(With $\overline{A}$ we denoted the closure of a set $A$). The hierarchy (7.19) is constructed such that $V$-spaces are self-similar,

$$f\left(2^j x\right) \in V_j\quad \mathrm{iff}\quad f(x) \in V_0{}.$$ (7.20)

with the requirement that there exists a scaling function $\phi \in V_0$ whose integer-translates span the space $V_0$,

$$V_0 = \left\{ f \in {{\mathbb{L}}}_2({{\mathbb{R}}})\vert ~~f(x)=\sum_k c_k \phi(x-k) \right\}{},$$

and for which the family $\{\phi(\bullet - k)$, $k \in {{\mathbb{Z}}}\}$ is an orthonormal basis. It can be assumed that $\int \phi(x) \mathrm{d} x \geq 0$. With this assumption this integral is in fact equal to $1$. Because of containment $V_0 \subset V_{1}$, the function $\phi(x) \in V_0$ can be represented as a linear combination of functions from $V_{1}$, i.e.,

$$\phi(x) = \sum_{k \in {{\mathbb{Z}}}} h_k \sqrt{2}\phi(2 x - k){},$$ (7.21)

for some coefficients $h_k$, $k \in {{\mathbb{Z}}}$. This equation called the scaling equation (or two-scale equation) is fundamental in constructing, exploring, and utilizing wavelets.

Theorem 2   For the scaling function it holds

$$\int_{{{\mathbb{R}}}} \phi(x) d x = 1{},$$

or, equivalently,

$$\Phi(0) = 1{},$$

where $\Phi(\omega)$ is Fourier transform of $\phi$, $\int_{{{\mathbb{R}}}} \phi(x) \mathrm{e}^{- \mathrm{i} \omega x} \mathrm{d} x$.

The coefficients $h_n$ in (7.21) are important in efficient application of wavelet transforms. The (possibly infinite) vector ${\boldsymbol{h}} = \{ h_n$, $n \in {{\mathbb{Z}}} \}$ will be called a wavelet filter. It is a low-pass (averaging) filter as will become clear later by its analysis in the Fourier domain.

To further explore properties of multiresolution analysis subspaces and their bases, we will often work in the Fourier domain.

It will be convenient to use Fourier domain for subsequent analysis of wavelet paradigm. Define the function $m_0$ as follows:

$$m_0(\omega) = \frac{1}{\sqrt{2}} \sum_{k \in {{\mathbb{Z}}}} h_k \mathrm{e}^{-\mathrm{i} k \omega} = \frac{1}{\sqrt{2}} H(\omega){}.$$ (7.22)

The function in (7.22) is sometimes called the transfer function and it describes the behavior of the associated filter $\boldsymbol{h}$ in the Fourier domain. Notice that the function $m_0$ is $2 \pi$-periodic and that filter taps $\{ h_n$, $n \in {{\mathbb{Z}}} \}$ are in fact the Fourier coefficients in the Fourier serias of $H(\omega)=\sqrt{2} m_0(\omega)$.

In the Fourier domain the relation (7.21) becomes

$$\Phi(\omega) = m_0 \left(\frac{\omega}{2}\right) \Phi\left(\frac{\omega}{2}\right){},$$ (7.23)

where $\Phi(\omega)$ is the Fourier transform of $\phi(x)$. Indeed,

$$\Phi(\omega) = \int_{-\infty}^{\infty} \phi(x) \mathrm{e}^{-\mathrm{i} \omega x} \mathrm{d} x \nonumber$$

$$= \sum_k \sqrt{2} h_k \int_{-\infty}^{\infty} \phi(2 x - k) \mathrm{e}^{-\mathrm{i} \omega x} \mathrm{d} x \nonumber$$

$$= \sum_k \frac{h_k}{\sqrt{2}} \mathrm{e}^{-\mathrm{i} k \omega/2}... ...2 x - k) \mathrm{e}^{-\mathrm{i} (2 x-k) \omega/2 } \mathrm{d}(2 x-k) \nonumber$$

$$= \sum_k \frac{ ~h_k}{\sqrt{2}} \mathrm{e}^{-\mathrm{i} k \omega/2}\Phi\left(\frac{\omega}{2}\right)\nonumber$$

$$= m_0\left(\frac{\omega}{2}\right) \Phi\left(\frac{\omega}{2}\right){}.$$

By iterating (7.23), one gets

$$\Phi(\omega) = \prod_{n=1}^{\infty} m_0\left(\frac{\omega}{2^n}\right){},$$ (7.24)

which is convergent under very mild conditions concerning the rates of decay of the scaling function $\phi$.

Next, we prove two important properties of wavelet filters associated with an orthogonal multiresolution analysis, normalization and orthogonality.

Normalization.

$$\index{normalization property}\sum_{k \in {{\mathbb{Z}}}} h_k = \sqrt{2}{}.$$ (7.25)

Proof:

$$\int \phi(x) \mathrm{d} x = \sqrt{2} \sum_k h_k \int \phi(2 x -k) \mathrm{d} x$$

$$= \sqrt{2} \sum_k h_k \frac{1}{2} \int \phi(2 x -k) d(2 x -k)$$

$$= \frac{\sqrt{2}}{2} \sum_k h_k \int \phi(x) \mathrm{d} x{}.$$

Since $\int \phi(x) d x \ne 0$ by assumption, (7.25) follows.

This result also follows from $m_0(0) = 1$.

Orthogonality.

For any $l \in {{\mathbb{Z}}},$

$$\index{orthogonality property}\sum_k h_{k} h_{k - 2 l} = \delta_{l}{}.$$ (7.26)

Proof: Notice first that from the scaling equation (7.21) it follows that

$$\phi(x) \phi(x-l) = \sqrt{2} \sum_k h_k \phi(2 x -k) \phi(x -l)$$ (7.27)

$$= \sqrt{2} \sum_k h_k \phi(2 x -k) \sqrt{2} \sum_m h_m \phi(2(x -l)-m){}.$$

By integrating the both sides in (7.27) we obtain

$$\delta_{l} = 2 \sum_k h_k \left[\sum_m h_m \frac{1}{2} \int \phi(2 x - k) \phi(2 x- 2 l - m) \, {\mathrm{d}}( 2 x ) \right]$$

$$= \sum_{k} \sum_m~ h_k h_m \delta_{k, 2 l +m}$$

$$= \sum_k h_k h_{k-2l}{}.$$

The last line is obtained by taking $k=2 l + m$.

An important special case is $l=0$ for which (7.26) becomes

$$\sum_k h^2_k = 1{}.$$ (7.28)

The fact that the system $\{\phi(\bullet - k)$, $k \in {{\mathbb{Z}}}\}$ constitutes an orthonormal basis for $V_0$ can be expressed in the Fourier domain in terms of either $\Phi(\omega)$ or  $m_0(\omega)$.

In terms of $\Phi(\omega)$,

$$\sum_{l=-\infty}^{\infty} \vert\Phi(\omega + 2 \pi l)\vert^2 = 1{}.$$ (7.29)

From the Plancherel identity and the $2 \pi$-periodicity of $\mathrm{e}^{\mathrm{i} \omega k}$ it follows

$$\delta_{k} = \int_{{{\mathbb{R}}}} \phi(x) \overline{ \phi(x-k)} \mathrm{d} x$$

$$= \frac{1}{2 \pi} \int_{{{\mathbb{R}}}} \Phi(\omega) \overline {\Phi(\omega)} \mathrm{e}^{\mathrm{i} \omega k} \mathrm{d} \omega$$

$$= \frac{1}{2 \pi} \int_0^{2 \pi} \sum_{l=-\infty}^{\infty} \vert\Phi(\omega+2 \pi l)\vert^2 \mathrm{e}^{\mathrm{i} \omega k} \mathrm{d} \omega{}.$$ (7.30)

The last line in (7.30) is the Fourier coefficient $a_k$ in the Fourier series decomposition of

$$f(\omega)= \sum_{l=-\infty}^{\infty} \vert\Phi(\omega+2 \pi l)\vert^2{}.$$

Due to the uniqueness of Fourier representation, $f(\omega)=1$. As a side result, we obtain that $\Phi(2 \pi n)= 0$, $n \ne 0$, and $\sum_n \phi(x - n) = 1$. The last result follows from inspection of coefficients $c_k$ in the Fourier decomposition of $\sum_n \phi(x - n)$, the series $\sum_k c_k \mathrm{e}^{2 \pi \mathrm{i} k x}$. As this function is 1-periodic,

$$c_k = \int_{0}^1 \left( \sum_n \phi(x-n) \right) \mathrm{e}^{- 2 ... ...athrm{e}^{- 2 \pi \mathrm{i} k x} \mathrm{d} x =\Phi(2 \pi k) = \delta_{0,k}{}.$$

Remark 1   Utilizing the identity (7.29), any set of independent functions spanning $V_0$, $\{ \phi(x - k)$, $k \in {{\mathbb{Z}}}\}$, can be orthogonalized in the Fourier domain. The orthonormal basis is generated by integer-shifts of the function

$$\mathcal{F}^{-1} \left[\frac{\Phi(\omega)}{\sqrt{ \sum_{l=-\infty}^{\infty} \vert\Phi(\omega+2 \pi l)\vert^2} }\right]{}.$$ (7.31)

This normalization in the Fourier domain is used in constructing of some wavelet bases.

Orthogonality condition 7.29 can be expressed in terms of $m_0$ as:

$$\left\vert m_0 (\omega) \right\vert^2 + \left\vert m_0 (\omega + \pi) \right\vert^2 = 1{}.$$ (7.32)

Since $\sum_{l=-\infty}^{\infty} \vert\Phi(2 \omega + 2 l \pi )\vert^2 =1$, then by (7.23)

$$\sum_{l=-\infty}^{\infty} \left\vert m_0(\omega + l \pi)\right\vert^2 \left\vert\Phi( \omega + l \pi )\right\vert^2 =1{}.$$ (7.33)

Now split the sum in (7.33) into two sums - one with odd and the other with even indices, i.e.,

$$1 = \sum_{k=-\infty}^{\infty} \left\vert m_0(\omega + 2 k \pi)\right\vert^2 \left\vert\Phi( \omega + 2 k \pi )\right\vert^2 +$$

$$\sum_{k=-\infty}^{\infty} \left\vert m_0(\omega + (2 k + 1) \pi)\right\vert^2 \left\vert\Phi( \omega + (2 k + 1) \pi )\right\vert^2{}.$$

To simplify the above expression, we use (7.29) and the $2 \pi$-periodicity of  $m_0(\omega)$.

$$1 = \vert m_0(\omega)\vert^2 \sum_{k=-\infty}^{\infty} \left\vert\P... ...2 \sum_{k=-\infty}^{\infty} \left\vert\Phi((\omega+\pi) + 2 k \pi)\right\vert^2$$

$$= \vert m_0(\omega)\vert^2 + \vert m_0(\omega + \pi)\vert^2{}.$$

Whenever a sequence of subspaces satisfies MRA properties, there exists (though not unique) an orthonormal basis for  ${{\mathbb{L}}}_2({{\mathbb{R}}})$,

$$\left\{ \psi_{jk}(x) = 2^{j/2} \psi\left( 2^{j} x - k\right){}, \quad j,k \in {{\mathbb{Z}}} \right\}$$ (7.34)

such that $\{\psi_{jk}(x)$, $j$-fixed, $k \in {{\mathbb{Z}}}\}$ is an orthonormal basis of the ''difference space'' $W_j = V_{j+1} \ominus V_j$. The function $\psi(x) = \psi_{00}(x)$ is called a wavelet function or informally the mother wavelet.

Next, we discuss the derivation of a wavelet function from the scaling function. Since $\psi(x) \in V_{1}$ (because of the containment $W_0 \subset V_{1}$), it can be represented as

$$\psi(x) = \sum_{k \in {{\mathbb{Z}}}} g_k \sqrt{2} \phi(2 x -k){},$$ (7.35)

for some coefficients $g_k, ~k \in {{\mathbb{Z}}}$.

Define

$$m_1(\omega) = \frac{1}{\sqrt{2}}~ \sum_k g_k \mathrm{e}^{-\mathrm{i} k \omega}{}.$$ (7.36)

By mimicking what was done with $m_0$, we obtain the Fourier counterpart of (7.35),

$$\Psi(\omega) = m_1\left(\frac{\omega}{2}\right) \Phi\left(\frac{\omega}{2}\right){}.$$ (7.37)

The spaces $W_0$ and $V_0$ are orthogonal by construction. Therefore,

$$0 = \int \psi(x) \phi(x -k) \mathrm{d} x = \frac{1}{2 \pi} \int \Psi(\omega) \overline { \Phi(\omega)} \mathrm{e}^{\mathrm{i} \omega k} \mathrm{d}\omega$$

$$= \frac{1}{2 \pi} \int_0^{2 \pi} \sum_{l=-\infty}^{\infty} \Psi(\... ...ne { \Phi(\omega+2 l \pi)} \mathrm{e}^{\mathrm{i} \omega k} \mathrm{d}\omega{}.$$

By repeating the Fourier series argument, as in (7.29), we conclude

$$\sum_{l=-\infty}^{\infty} \Psi(\omega+2 l \pi) \overline { \Phi(\omega+2 l \pi)} = 0{}.$$

By taking into account the definitions of $m_0$ and $m_1$, and by the derivation as in (7.32), we find

$$m_1(\omega) \overline{ m_0(\omega) } + m_1(\omega+\pi) \overline{ m_0(\omega+\pi) } = 0{}.$$ (7.38)

From (7.38), we conclude that there exists a function  $\lambda(\omega)$ such that

$$\left(m_1(\omega), m_1(\omega+\pi)\right) = \lambda(\omega) \left( \overline{ m_0(\omega+\pi) }, - \overline{ m_0(\omega)} \right){}.$$ (7.39)

By substituting $\xi = \omega + \pi$ and by using the $2 \pi$-periodicity of $m_0$ and $m_1$, we conclude that

$$\lambda (\omega) = - \lambda(\omega + \pi){},\quad \mathrm{and}$$ (7.40)

$$\lambda (\omega)\quad \mathrm{is} 2 \pi\mathrm{-periodic}{}.$$

Any function $\lambda(\omega)$ of the form $\mathrm{e}^{\pm \mathrm{i} \omega} S(2 \omega)$, where $S$ is an ${{\mathbb{L}}}_2([0,2\pi])$, $2 \pi$-periodic function, will satisfy (7.38); however, only the functions for which $\vert\lambda(\omega)\vert=1$ will define an orthogonal basis $\psi_{jk}$ of  ${{\mathbb{L}}}_2({{\mathbb{R}}})$.

To summarize, we choose $\lambda(\omega)$ such that

(i)

$\lambda(\omega)$ is $2 \pi$-periodic,

(ii)

$\lambda(\omega) = -\lambda(\omega + \pi)$, and

(iii)

$\vert\lambda(\omega)\vert^2 = 1$.

Standard choices for $\lambda(\omega)$ are $-\mathrm{e}^{-\mathrm{i} \omega}$, $\mathrm{e}^{-\mathrm{i} \omega}$, and $\mathrm{e}^{\mathrm{i} \omega}$; however, any other function satisfying (i)-(iii) will generate a valid $m_1$. We choose to define $m_1(\omega)$ as

$$m_1(\omega) = - \mathrm{e}^{-\mathrm{i}\omega} \overline{ m_0(\omega + \pi)}{}.$$ (7.41)

since it leads to a convenient and standard connection between the filters  $\boldsymbol{h}$ and  $\boldsymbol{g}$.

The form of $m_1$ and (7.29) imply that $\{\psi(\bullet - k)$, $k \in {{\mathbb{Z}}}\}$ is an orthonormal basis for $W_0$.

Since $\vert m_1(\omega)\vert=\vert m_0(\omega+\pi)\vert$, the orthogonality condition (7.32) can be rewritten as

$$\vert m_0(\omega)\vert^2+\vert m_1(\omega)\vert^2 = 1{}.$$ (7.42)

By comparing the definition of $m_1$ in (7.36) with

$$m_1(\omega) =-\mathrm{e}^{-\mathrm{i} \omega} \frac{1}{\sqrt{2}} \sum_k h_k \mathrm{e}^{\mathrm{i} (\omega + \pi) k}$$

$$= \frac{1}{\sqrt{2}} \sum_k (-1)^{1-k} h_k \mathrm{e}^{-\mathrm{i} \omega (1-k)}$$

$$= \frac{1}{\sqrt{2}} \sum_n (-1)^n h_{1-n} \mathrm{e}^{-\mathrm{i} \omega n}{},$$

we relate $g_n$ and $h_n$ as

$$g_n = (-1)^n ~h_{1-n}{}.$$ (7.43)

In signal processing literature, (7.43) is known as the quadrature mirror relation and the filters  $\boldsymbol{h}$ and  $\boldsymbol{g}$ as quadrature mirror filters.

Remark 2   Choosing $\lambda(\omega)=\mathrm{e}^{\mathrm{i}\omega}$ leads to the rarely used high-pass filter $g_n = (-1)^{n-1} h_{-1-n}$. It is convenient to define $g_n$ as $(-1)^n h_{1-n + M}$, where $M$ is a ''shift constant.'' Such re-indexing of  $\boldsymbol{g}$ affects only the shift-location of the wavelet function.

7.3.4 Haar Wavelets

In addition to their simplicity and formidable applicability, Haar wavelets have tremendous educational value. Here we illustrate some of the relations discussed in the Sect. 7.3.3 using the Haar wavelet. We start with scaling function $\phi(x)=\boldsymbol{1}(0 \leq x \leq 1)$ and pretend that everything else is unknown. By inspection of simple graphs of two scaled Haar wavelets $\phi(2 x)$ and $\phi(2 x + 1)$ stuck to each other, we conclude that the scaling equation (7.21) is

$$\phi(x) =\phi(2 x)+ \phi(2 x -1)$$

$$=\frac{1}{\sqrt{2}} \sqrt{2} \phi(2 x) +\frac{1}{\sqrt{2}} \sqrt{2} \phi(2 x-1){},{}$$ (7.44)

which yields the wavelet filter coefficients:

$$h_0 = h_1= \frac{1}{\sqrt{2}}{}.$$

The transfer functions are

$$m_0(\omega) = \frac{1}{\sqrt{2}} \left(\frac{1}{\sqrt{2}} \mathrm... ...\mathrm{i} \omega 1}\right) = \frac{1 + \mathrm{e}^{- \mathrm{i} \omega}}{2}{}.$$

and

$$m_1(\omega) = - \mathrm{e}^{-\mathrm{i} \omega} \overline{m_0 (\o... ...}^{ \mathrm{i} \omega}\right)= \frac{1 - \mathrm{e}^{-\mathrm{i} \omega}}{2}{}.$$

Notice that $m_0(\omega) = \vert m_0(\omega)\vert \mathrm{e}^{\mathrm{i} \varphi(\omega)} = \cos(\omega/2) \cdot \mathrm{e}^{-\mathrm{i} \omega/2}$ (after $\cos x = (\mathrm{e}^{\mathrm{i} x}+\mathrm{e}^{-\mathrm{i} x})/2$). Since $\varphi(\omega) = - \frac{\omega}{2}$, the Haar wavelet has linear phase, i.e., the scaling function is symmetric in the time domain. The orthogonality condition $\vert m_0(\omega)\vert^2 + \vert m_1(\omega)\vert^2 = 1$ is easily verified, as well.

Relation (7.37) becomes

$$\Psi(\omega) = \frac{1 - \mathrm{e}^{-\mathrm{i} \omega/2}}{2} \P... ...ac{1}{2} \Phi\left(\frac{\omega}{2}\right) \mathrm{e}^{-\mathrm{i} \omega/2}{},$$

and by applying the inverse Fourier transform we obtain

$$\psi(x) = \phi(2 x) - \phi(2 x -1)$$

in the time-domain. Therefore we ''have found'' the Haar wavelet function $\psi$. From the expression for $m_1$ or by inspecting the representation of $\psi(x)$ by $\phi(2 x)$ and $\phi(2 x -1)$, we ''conclude'' that $g_0=-g_{-1} = \frac{1}{\sqrt{2}}$.

Although the Haar wavelets are well localized in the time domain, in the frequency domain they decay at the slow rate of $O(1/n)$ and are not effective in approximating smooth functions.

7.3.5 Daubechies' Wavelets

The most important family of wavelets was discovered by Ingrid Daubechies and fully described in Daubechies (1992)[8]. This family is compactly supported with various degrees of smoothness.

The formal derivation of Daubechies' wavelets goes beyond the scope of this chapter, but the filter coefficients of some of its family members can be found by following considerations.

For example, to derive the filter taps of a wavelet with $N$ vanishing moments, or equivalently, $2 N$ filter taps, we use the following equations.

The normalization property of scaling function implies

$$\sum_{i = 0}^{2 N -1} h_i = \sqrt{2}{},$$

requirement for vanishing moments for wavelet function $\psi$ leads to

$$\sum_{i = 0}^{2 N -1} (-1)^i i^k h_i = 0{}, \quad k=0,1,\ldots,N-1{},$$

and, finally, the orthogonality property can be expressed as

$$\index{orthogonality property} \sum_{i = 0}^{2 N -1} h_i h_{i+2k} = \delta_k\quad k=0,1,\ldots,N-1{}.$$

We obtained $2 N+ 1$ equations with $2 N$ unknowns; however the system is solvable since the equations are not linearly independent.

Example 4   For $N=2$, we obtain the system:

$$\left\{ \begin{array}{l} \ h_0 + h_1 + h_2 + h_3 = \sqrt{2} \\ \ ... ... - h_1 + 2 h_2 - 3 h_3 = 0 , \\ \ h_0~ h_2 + h_1~h_3 = 0 \end{array}{}, \right.$$

which has a solution $h_0= \frac{1+\sqrt{3}}{4 \sqrt{2}}$, $h_1= \frac{3+\sqrt{3}}{4 \sqrt{2}}$, $h_2= \frac{3-\sqrt{3}}{4 \sqrt{2}}$, and $h_3= \frac{1-\sqrt{3}}{4 \sqrt{2}}$.

For $N=4$, the system is

$$\left\{ \begin{array}{l} \ h_0 + h_1 + h_2 + h_3 + h_4 + h_5 + h_... ...h_2 - 27 h_3 + 64 h_4 - 125 h_5 + 216 h_6 - 343 h_7 = 0 {}. \end{array} \right.$$

Figure 7.9 depicts two scaling function and wavelet pairs from the Daubechies family. Figure 7.9a,b depict the pair with two vanishing moments, while Fig. 7.9c,d depict the pair with four vanishing moments.

Figure 7.9: Wavelet functions from Daubechies family. (a) Daubechies scaling function, $2$ vanishing moments, $4$ tap filter (b) Wavelet function corresponding to (a), (c) Daubechies scaling function, $4$ vanishing moments, $8$ tap filter (d) Wavelet function corresponding to (c)

(a)

(b)

(c)

(d)