Non-decimated Wavelet Transform for a Shift-invariant Analysis

ABSTRACT

Due to the ability of time-frequency location, the wavelet transform has been applied in several areas of research involving signal analysis and processing, often replacing the conventional Fourier transform. The discrete wavelet transform has great application potential, being an important tool in signal compression, signal and image processing, smoothing and de-noising data. It also presents advantages over the continuous version because of its easy implementation, good computational performance and perfect reconstruction of the signal upon inversion. Nevertheless, the downsampling required in the computation of the discrete wavelet transform makes it shift variant and not appropriated to some applications, such as for signals or time series analysis. On the other hand, the Non-Decimated Discrete Wavelet Transform is shift-invariant because it eliminates the downsampling and, consequently, is more appropriate for identifying both stationary and non-stationary behaviors in signals. However, the non-decimated wavelet transform has been underused in the literature. This paper intends to show the advantages of using the non-decimated wavelet transform in signal analysis. The main theoretical and practical aspects of the multi-scale analysis of time series from non-decimated wavelets in terms of its formulation using the same pyramidal algorithm of the decimated wavelet transform was presented. Finally, applications with a simulated and real time series compare the performance of the decimated and non-decimated wavelet transform, demonstrating the superiority of non-decimated one, mainly due to the shift-invariant analysis, patterns detection and more perfect reconstruction of a signal.

Keywords:
Non-decimated wavelets; shift invariance; time series; signal analysis

RESUMO

Devido sua habilidade de localização tempo-frequência, a Transformada wavelet tem sido aplicada em várias áreas de pesquisa envolvendo análise e processamento de dados, frequentemente substituindo a convencional Transformada de Fourier. A TransformadaWavelet Discreta tem um grande potencial de aplicação¸ destacando-se como uma importante ferramenta na compressão de sinal, processamento de imagem e sinal, suavização e filtragem de ruídos em dados. Ela também apresenta vantagens sobre a versão contínua por causa de sua fácil implementação, bom desempenho computacional e reconstrução perfeita do sinal após inversão. No entanto, a decimação requerida no cálculo da Transformada Wavelet Discreta a torna variante à translação e não apropriada para algumas aplicaões, tais como análise de sinais ou séries temporais. Por outro lado, a Transformada Wavelet Discreta Não Decimada é invariante à translação, porque elimina o processo de decimação, e consequentemente, é mais apropriada para identificar comportamentos estacionários e não estacionários presentes no sinal. No entanto, a Transformada Wavelet Não Decimada tem sido pouco usada na literatura. Esse artigo pretende mostrar as vantagens do uso na Transformada Wavelet Não Decimada na análise de sinais. Os principais aspectos teóricos e práticos da análise multiescala de séries temporais a partir das wavelets não decimadas, em termos de sua formulação usando o mesmo algoritmo piramidal da Transformada Wavelet Decimada, são apresentados. Por fim, aplicações com séries temporais simuladas e reais comparam o desempenho das transformadas wavelet decimada e não decimada, demonstrando a superioridade da wavelet não decimada, principalmente devido à análise invariante a translação, detecção de padrões e uma reconstrução mais perfeita do sinal.

Palavras-chave:
Wavelets não decimadas; invariância à translação; séries temporais; análise de sinais

1 INTRODUCTION

The wavelet transform allows extracting information of stationary and non-stationary signal variations in time and frequency, i.e., identifying their frequency of occurrence, localization in time, and making a reliable approximation of magnitude of this variation. For this reason, this methodology has been adopted for a vast number of applications, often replacing the conventional Fourier transform.

There are two kinds of wavelet transforms, the continuous wavelet transform (CWT) ¹¹ and the discrete wavelet transforms, with its decimated (DWT) ⁶ and non-decimated (NDWT) ⁷ versions. Both wavelet transforms differ in the representation of the scale and location parameters of the wavelet function, which can take continuous or discrete (powers of two) values, respectively.

These differences result in advantages and disadvantages for the two classes of wavelet transforms and also determine in which applications each wavelet transform can provide superior results.

The DWT, for example, provides a sparse representation for many natural signals, it is therefore an important tool in signal compression. The DWT is an orthonormal transform (when using an orthogonal wavelet) able to separate the signal from the noise. This is possible because the signal of interest is typically captured by a few large-magnitude DWT coefficients, while the noise results in many small DWT coefficients, which can be throw away without harming the quality of the signal approximation ⁷. As a result, the important features of the signal are captured by a subset of DWT coefficients that is typically much smaller than the original signal, namely, the compressed signal. The same considerations are taken in noise filtering or de-noising. The CWT, instead, is a highly redundant transform and it not appropriate to these applications. The computational resources required to compute the CWT and store the coefficients is much larger than the DWT.

On the other hand, because the DWT downsamples, a shift in the input signal does not manifest itself as an equivalent shift in the DWT coefficients at all levels, i.e, the DWT is not shift-invariant. Thus, a simple shift in a signal can cause a significant change of signal energy in the DWT coefficients by scale, what is not ideal for some applications. The CWT and NDWT, instead, are shift-invariant so they are perfect to time series analysis.

The other advantage of the discrete wavelet transforms in relation to CWT is the easy implementation and good computational performance. The discrete approach of the wavelet transform can be performed with Mallat’s and the ‘à trous’ algorithms ⁹. The first is an orthogonal, dyadic, non-symmetric, decimated and non-redundant algorithm. The ’à trous’, in opposite, is a non-orthogonal, dyadic, symmetric, shift-invariant and redundant algorithm ⁴. Both algorithms are equivalent to discrete filter banks, where the signal is iteratively filtered by a low-pass and a high-pass filter. The Mallat’s algorithm, also called the pyramidal algorithm, is the most used to computation of the discrete wavelet transform. In the DWT, the filter outputs are downsampled at each successive stage of the pyramidal algorithm, namely, for each two outputs of the filter, one output is discarded. In the NDWT, however, the outputs are not downsampled, wherein each scale will have the same number of wavelet coefficients. The filters that define the discrete wavelet transforms typically only have a small number of coefficients so the transform can be implemented very efficiently. Furthermore, for the most common implementation of the CWT it is necessary that the wavelet is defined in a closed-form expression, while for both DWT and NDWT, only the filters are sufficient.

Another key point of the discrete wavelet transform is the perfect reconstruction of the signal upon inversion. An inverse CWT can be implemented, but often the reconstructions are not perfect. Reconstructing a signal from the CWT coefficients is a much less stable numerical operation.

It should be emphasized, however, that the finer sampling of scales in the CWT typically results in a higher-fidelity signal analysis. Continuous analysis is often easier to interpret because its redundancy tends to reinforce the traits and makes all information more visible. Thus, it is possible to localize transients in a signal or characterize periodic behavior better with the CWT than with the discrete wavelet transforms. Nevertheless, it is possible to have a discrete wavelet transform, which does not lose important information and has the advantages of implementation and computational effort. This is the case of the NDWT, which can be seen as a compromise between the DWT and CWT because of its redundancy, but not as redundant as CWT. Although, neither the CWT nor the NDWT are orthonormal transforms, the NDWT can be computed similarly to the ordinary DWT but without downsampling, ensuring the shift invariance, which is ideal for analyzing time series.

While it is well known that the DWT is efficient to data compression, signal and image processing, smoothing and de-noising data, among others applications; for signals or time series analysis, however, the NDWT shows to be more appropriate. As the NDWT is shift-invariant and represents a time series with the same number of coefficients at each scale, it is possible to detect the occurrence of hidden information such as stationary or non-stationary patterns and its time/frequency location. Recently, the NDWT has been applied in some areas of research: in GNSS signal analysis to investigate the ionospheric scintillation effect ²; in daily temperature data to explore the time scale patterns in the relationship between the average ambient temperature and the number of deaths due to cardiovascular diseases ¹; in analysis of BSE and NSE indexes financial time series ⁵; in waves pressure, height geopotential and thickness time series analysis in order to forecast rain events ¹⁰, among other applications. However, the use of NDWT should be much larger in signal analysis. Thus, in this paper we aim to point out the advantages of NDWT for signal or time series analysis in an intuitive point of view of Mallat´s pyramidal algorithm. As this algorithm is well known and largely used in the literature, the implementation of NDWT from it can be facilitated. The focus in this paper is in the sense of evaluation of methodologies that allow the signal analysis with possibility of investigating behavior in it as well as the perfect reconstruction. Thus, only discrete transforms will be taken into account, i.e. both decimated (DWT) and non-decimated (NDWT) transforms. Simulated and real signals are used to illustrate the results.

2 NDWT SIGNAL ANALYSIS

2.1 NDWT wavelet and scaling filters

The NDWT of a time series X = [X ₀, X ₁,...,X _n-1 ] in a level J yields column vectors ${\tilde{W}}_1,{\tilde{W}}_2\mathrm{,...,}{\tilde{W}}_J$ and ${\tilde{V}}_J$ of length n. The vector ${\tilde{W}}_j$ contains the NDWT wavelet coefficients associated with changes in X on a scale of τ _j = 2^j-1 , j = 1,...,J, while ${\tilde{V}}_J$ contains the NDWT scaling coefficients associated with variations at scale τ _j = 2 ^{J (8}. Such wavelet and scaling coefficients are the filtering result on the time series X with the NDWT wavelet and scaling filters.

A wavelet filter {h _k } must satisfy the following properties:

1) ${\displaystyle {\sum }_{k=0}^{K-1}{h}_k=0}$;

2) ${\displaystyle {\sum }_{k=0}^{K-1}{h}_k^2=1}$;

3) ${\displaystyle {\sum }_{k=0}^{K-1}{h}_k{h}_{k+2n}={\displaystyle {\sum }_{-\infty }^{+\infty }{h}_k{h}_{k+2n}=0}}$, for most of nonzero integer K values and all nonzero integers n.

On other words, a wavelet filter must sum zero, have unit energy, and must be orthogonal to its even shifts.

The second required filter {g _k } is obtained by quadrature mirror filter $g_k={(-1)}^{k+1}{h}_{K-1-k}$. The filter {g _k } is known as the scaling filter and must also satisfy three basic properties:

1) ${\displaystyle {\sum }_{k=0}^{K-1}{g}_k=\sqrt{2}}$;

2) ${\displaystyle {\sum }_{k=0}^{K-1}{g}_k^2=1}$;

3) ${\displaystyle {\sum }_{k=0}^{K-1}{g}_k{g}_{k+2n}={\displaystyle {\sum }_{-\infty }^{+\infty }{h}_k{h}_{k+2n}=0}}$, for all nonzero integers n.

The NDWT wavelet filter $\left\{{\tilde{h}}_k\right\}$ and NDWT scaling filter $\left\{{\tilde{g}}_k\right\}$ are the rescaled version of the wavelet filter {h _k } and scaling {g _k } shown above, defined via ${\tilde{h}}_k=h_k/\sqrt{2}$ and ${\tilde{g}}_k=g_k/\sqrt{2}$.

So, the result of the filtering of a time series {X _t: t = 0, …, N - 1} with the NDWT wavelet and scaling filters is given, respectively, by

$\begin{array}{l}{\tilde{W}}_{\mathrm{1,}t}={\displaystyle {\sum }_{k=0}^{K-1}{\tilde{h}}_k{X}_{t-k\mathrm{mod}n}}\\ {\tilde{V}}_{\mathrm{1,}t}={\displaystyle {\sum }_{k=0}^{K-1}{\tilde{g}}_k{X}_{t-k\mathrm{mod}n},t=\mathrm{0,1,...,}n-1.}\end{array}$ (2.1)

These two sequences are the NDWT in the level J = 1. The term mod in (2.1) allows a circularly filtering, consequently, X is represented with the same number of coefficients at each scale.

A relationship between the DWT and NDWT wavelet and scaling coefficients, can be expressed

$\begin{array}{l}{W}_{\mathrm{1,}t}\equiv 2^{1/2}{\tilde{W}}_{\mathrm{1,2}t+1}={\displaystyle {\sum }_{k=0}^{K-1}{\tilde{h}}_k{X}_{2t+1-k\mathrm{mod}n}}\\ V_{\mathrm{1,}t}\equiv 2^{1/2}{\tilde{V}}_{\mathrm{1,2}t+1}={\displaystyle {\sum }_{k=0}^{K-1}{\tilde{g}}_k{X}_{2t+1-k\mathrm{mod}n},}t=\mathrm{0,...,}\frac{n}{2}-1.\end{array}$ (2.2)

2.2 NDWT Formulation

The great motivation for formulating the NDWT is to define a transform that acts similarly to the DWT, but does not suffer sensitivity to the choice of a starting point for a time series, in other words, it is shift-invariant.This DWT’sensitivity is entirely due to downsampling the outputs from the wavelet and scaling filters at each stage of the pyramidal algorithm. To make the understanding easier, we are going to use a matrix notation ⁸. The filter outputs that usually discarded at first stage of the DWT pyramidal algorithm can be obtained by applying the DWT pyramidal algorithm to the circularly shifted vector TX rather than X. It is worth mentioning that if the vector X=[X ₀,X ₁,…,X _n-1 ]^t , so TX=[X _n-1 ,X ₀,…,X _n-2 ]^t is the vector X shifted by one unit, where T _n is called translation matrix, given by

$T=\left[\begin{array}{cccccccc}0& 0& 0& 0& \cdots & 0& 0& 1\\ 1& 0& 0& 0& \cdots & 0& 0& 0\\ 0& 1& 0& 0& \cdots & 0& 0& 0\\ ⋮& ⋮& ⋮& ⋮& \cdots & ⋮& ⋮& ⋮\\ 0& 0& 0& 0& \cdots & 1& 0& 0\\ 0& 0& 0& 0& \cdots & 0& 1& 0\end{array}\right].$ (2.3)

This procedure suggests how to eliminate the downsampling and defines the first stage of NDWT pyramidal algorithm when n is an even sample size. The idea is to apply the usual DWT pyramidal algorithm twice, once to X and once to TX, and then merging the two sets of DWT coefficients together ⁸. The first application yields in

$\left[\begin{array}{c}{W}_1\\ V_1\end{array}\right]=\left[\begin{array}{c}{B}_1\\ A_1\end{array}\right]X=P_{1}X,$ (2.4)

where,

$B_1=\left[\begin{array}{ccccccccccccc}{h}_1& h_0& 0& 0& 0& \cdots & 0& 0& 0& 0& 0& h_3& h_2\\ h_3& h_2& h_1& h_0& 0& \cdots & 0& 0& 0& 0& 0& 0& 0\\ ⋮& ⋮& ⋮& ⋮& ⋮& & ⋮& ⋮& ⋮& ⋮& ⋮& ⋮& ⋮\\ 0& 0& 0& 0& 0& \cdots & 0& h_3& h_2& h_1& h_0& 0& 0\\ 0& 0& 0& 0& 0& \cdots & 0& 0& 0& h_3& h_2& h_1& h_0\end{array}\right],$

$A_1=\left[\begin{array}{ccccccccccccc}{g}_1& g_0& 0& 0& 0& \cdots & 0& 0& 0& 0& 0& g_3& g_2\\ g_3& g_2& g_1& g_0& 0& \cdots & 0& 0& 0& 0& 0& 0& 0\\ ⋮& ⋮& ⋮& ⋮& ⋮& & ⋮& ⋮& ⋮& ⋮& ⋮& ⋮& ⋮\\ 0& 0& 0& 0& 0& \cdots & 0& g_3& g_2& g_1& g_0& 0& 0\\ 0& 0& 0& 0& 0& \cdots & 0& 0& 0& g_3& g_2& g_1& g_0\end{array}\right]$

and

$P_1=\left[\begin{array}{c}{B}_1\\ A_{!}\end{array}\right].$

In displaying the elements of the matrices, we specialize to the case K=4 and N>K for clarity, but the mathematical treatment holds in general.

In view of equation (2.2), we can denote the elements of W ₁ and V ₁ by

$\begin{array}{l}{W}_1={\left[2^{1/2}{\tilde{W}}_{\mathrm{1,1}}{\mathrm{,2}}^{1/2}{\tilde{W}}_{\mathrm{1,3}}{\mathrm{,...,2}}^{1/2}{\tilde{W}}_{\mathrm{1,}n-1}\right]}^t\\ V_1={\left[2^{1/2}{\tilde{V}}_{\mathrm{1,1}}{\mathrm{,2}}^{1/2}{\tilde{V}}_{\mathrm{1,3}}{\mathrm{,...,2}}^{1/2}{\tilde{V}}_{\mathrm{1,}n-1}\right]}^t.\end{array}$

Note that W ₁ and V ₁ contain all the odd indexed elements of the length n sequences $\left\{2^{1/2}{\tilde{W}}_{\mathrm{1,}t}\right\}$ and $\left\{2^{1/2}{\tilde{V}}_{\mathrm{1,}t}\right\}$, which are formed by circularly convolving the time series X with, respectively, the wavelet filter {h _k } and scaling filter {g _k }.

The second application consists of replacing X for TX and apply the DWT to the circularly shifted vector. Therefore, it is obtained

$\left[\begin{array}{c}{W}_{T\mathrm{,1}}\\ V_{T\mathrm{,1}}\end{array}\right]=P_{1}TX.$

Defining

$P_{T\mathrm{,1}}=P_{1}T=\left[\begin{array}{c}{B}_1\\ A_1\end{array}\right]T=\left[\begin{array}{c}{B}_{1}T\\ A_{1}T\end{array}\right]=\left[\begin{array}{c}{B}_{T\mathrm{,1}}\\ A_{T\mathrm{,1}}\end{array}\right],$

we can write

$\left[\begin{array}{c}{W}_{T\mathrm{,1}}\\ V_{T\mathrm{,1}}\end{array}\right]=P_{T\mathrm{,1}}X=\left[\begin{array}{c}{B}_{T\mathrm{,1}}\\ A_{T\mathrm{,1}}\end{array}\right]X,$ (2.5)

where

$B_{T\mathrm{,1}}=\left[\begin{array}{ccccccccccccc}{h}_0& 0& 0& 0& 0& \cdots & 0& 0& 0& 0& h_3& h_2& h_1\\ h_2& h_1& h_0& 0& 0& \cdots & 0& 0& 0& 0& 0& 0& h_3\\ ⋮& ⋮& ⋮& ⋮& ⋮& & ⋮& ⋮& ⋮& ⋮& ⋮& ⋮& ⋮\\ 0& 0& 0& 0& 0& \cdots & h_3& h_2& h_1& h_0& 0& 0& 0\\ 0& 0& 0& 0& 0& \cdots & 0& 0& h_3& h_2& h_1& h_0& 0\end{array}\right]$

and A _T,1 has the same structure as the above with h _k replaced by g _k . As B ₁ X is formed by odd indexed values of the sequence $\left\{2^{1/2}{\tilde{W}}_{\mathrm{1,}t}\right\}$, then, B _T,1 X, which is B ₁ X one unit shifted, will be formed of the even indexed values of the filter output $\left\{2^{1/2}{\tilde{W}}_{\mathrm{1,}t}\right\}$, that is

$W_{T\mathrm{,1}}={\left[2^{1/2}{\tilde{W}}_{\mathrm{1,0}}{\mathrm{,2}}^{1/2}{\tilde{W}}_{\mathrm{1,2}}{\mathrm{,...,2}}^{1/2}{\tilde{W}}_{\mathrm{1,}n-2}\right]}^t,$

and by a similar argument the elements of V _T,1 are given by

$V_{T\mathrm{,1}}={\left[2^{1/2}{\tilde{V}}_{\mathrm{1,0}}{\mathrm{,2}}^{1/2}{\tilde{V}}_{\mathrm{1,2}}{\mathrm{,...,2}}^{1/2}{\tilde{V}}_{\mathrm{1,}n-2}\right]}^t.$

So, it is possible to form the NDWT wavelet coefficients ${\tilde{W}}_1$ by rescaling the interleaved elements of W ₁ and W _1,T , and similarly obtain NDWT scaling coefficients ${\tilde{V}}_1$ from V ₁ and V _T,1 , i.e.

$\begin{array}{l}{\tilde{W}}_1={\left[{\tilde{W}}_{\mathrm{1,0}},{\tilde{W}}_{\mathrm{1,1}},{\tilde{W}}_{\mathrm{1,2}}\mathrm{,...,}{\tilde{W}}_{\mathrm{1,}n-1}\right]}^t\\ {\tilde{V}}_1={\left[{\tilde{V}}_{\mathrm{1,0}},{\tilde{V}}_{\mathrm{1,1}},{\tilde{V}}_{\mathrm{1,2}}\mathrm{,...,}{\tilde{V}}_{\mathrm{1,}n-1}\right]}^t.\end{array}$ (2.6)

Note in (2.6), that the elements of ${\tilde{W}}_1$ and ${\tilde{V}}_1$ are exactly the filters outputs ${\tilde{W}}_{\mathrm{1,}t}$ and ${\tilde{V}}_{\mathrm{1,}t}$ obtained by using the NDWT filters ${\tilde{h}}_k$ and ${\tilde{g}}_k$, as can be seen in (2.1).

Defining ${\tilde{B}}_1$ as the N?N matrix formed by interleaving the rows of B _T,1 and B ₁ and replacing each h _k by ${\tilde{h}}_k$, i.e.,

${\tilde{B}}_1\equiv \left[\begin{array}{ccccccccccc}{\tilde{h}}_0& 0& 0& \cdots & 0& 0& 0& 0& {\tilde{h}}_3& {\tilde{h}}_2& {\tilde{h}}_1\\ {\tilde{h}}_1& {\tilde{h}}_0& 0& \cdots & 0& 0& 0& 0& 0& {\tilde{h}}_3& {\tilde{h}}_2\\ {\tilde{h}}_2& {\tilde{h}}_1& {\tilde{h}}_0& \cdots & 0& 0& 0& 0& 0& 0& {\tilde{h}}_3\\ ⋮& ⋮& ⋮& ⋮& ⋮& & ⋮& ⋮& ⋮& ⋮& ⋮\\ 0& 0& 0& \cdots & 0& {\tilde{h}}_3& {\tilde{h}}_2& {\tilde{h}}_1& {\tilde{h}}_0& 0& 0\\ 0& 0& 0& \cdots & 0& 0& {\tilde{h}}_3& {\tilde{h}}_2& {\tilde{h}}_1& {\tilde{h}}_0& 0\\ 0& 0& 0& \cdots & 0& 0& 0& {\tilde{h}}_3& {\tilde{h}}_2& {\tilde{h}}_1& {\tilde{h}}_0\end{array}\right],$

we obtain ${\tilde{W}}_1={\tilde{B}}_{1}X$. With an analogous definition for ${\tilde{A}}_1$, we have ${\tilde{V}}_1={\tilde{A}}_{1}X$. Lastly, it is possible to represent the first stage of the NDWT pyramid algorithm as

$\left[\begin{array}{c}{\tilde{W}}_1\\ {\tilde{V}}_1\end{array}\right]=\left[\begin{array}{c}{\tilde{B}}_1\\ {\tilde{A}}_1\end{array}\right]X={\tilde{P}}_{1}X,$

where

${\tilde{P}}_1=\left[\begin{array}{c}{\tilde{B}}_1\\ {\tilde{A}}_1\end{array}\right].$

Because $P_1^t{P}_1=l_n$ and T ^t T=l _n , where T is a translation matrix given in (2.3), it follows that

$P_{T\mathrm{,1}}^t{P}_{T\mathrm{,1}}=T^t{P}_1^t{P}_{1}T=l_n$

and therefore P _T,1 is an orthonormal matrix. Hence, we obtain the following decompositions to X

$\begin{array}{l}{\Vert X\Vert }^2={\Vert W_1\Vert }^2+{\Vert V_1\Vert }^2\\ {\Vert X\Vert }^2={\Vert W_{T\mathrm{,1}}\Vert }^2+{\Vert V_{T\mathrm{,1}}\Vert }^2.\end{array}$

Furthermore, since

$\begin{array}{l}{\Vert W_1\Vert }^2+{\Vert W_{T\mathrm{,1}}\Vert }^2=2{\Vert {\tilde{W}}_1\Vert }^2\\ {\Vert V_1\Vert }^2+{\Vert V_{T\mathrm{,1}}\Vert }^2=2{\Vert {\tilde{V}}_1\Vert }^2.\end{array}$

it follows that

${\Vert X\Vert }^2={\Vert {\tilde{W}}_1\Vert }^2+{\Vert {\tilde{V}}_1\Vert }^2.$

From equations (2.4) and (2.5) has

$X=\left[B_1^t,A_1^t\right]\left[\begin{array}{c}{W}_1\\ V_1\end{array}\right]$

and

$X=\left[B_{T\mathrm{,1}}^t,A_{T\mathrm{,1}}^t\right]\left[\begin{array}{c}{W}_{T\mathrm{,1}}\\ V_{T\mathrm{,1}}\end{array}\right].$

Thus, X can be written as

$\begin{array}{l}X=\frac{1}{2}\left[B_1^t,A_1^t\right]\left[\begin{array}{c}{W}_1\\ V_1\end{array}\right]+\frac{1}{2}\left[B_{T\mathrm{,1}}^t,A_{T\mathrm{,1}}^t\right]\left[\begin{array}{c}{W}_{T\mathrm{,1}}\\ V_{T\mathrm{,1}}\end{array}\right]\\ \begin{array}{c}\begin{array}{c}\end{array}\end{array}\begin{array}{c}\end{array}=\frac{1}{2}\left(B_1^t{W}_1+A_1^t{V}_1+B_{T\mathrm{,1}}^t{W}_{T\mathrm{,1}}+A_{t\mathrm{,1}}^t{V}_{T\mathrm{,1}}\right)\\ \begin{array}{c}\end{array}\begin{array}{c}\end{array}=\frac{1}{2}\left(B_1^t{W}_1+B_{T\mathrm{,1}}^t{W}_{T\mathrm{,1}}\right)+\frac{1}{2}\left(A_1^t{V}_1+A_{T\mathrm{,1}}^t{V}_{T\mathrm{,1}}\right)\\ \begin{array}{c}\end{array}\begin{array}{c}\end{array}={\tilde{B}}_1^t{\tilde{W}}_1+{\tilde{A}}_1^t{\tilde{V}}_1\\ \begin{array}{c}\end{array}\begin{array}{c}\end{array}={\tilde{D}}_1+{\tilde{S}}_1\end{array}$

where ${\tilde{D}}_1\equiv {\tilde{B}}_1^t{\tilde{W}}_1$ is the first level NDWT detail and ${\tilde{S}}_1{\tilde{A}}_1^t{\tilde{V}}_1$ is the smooth corresponding.

We obtained until now the first level of NDWT coefficients. In the next sections, such concepts will be generalized in order to define the level j of the NDWT. It will also be presented the NDWT pyramidal algorithm, in order to calculate the wavelet and scaling coefficients in j levels of details.

2.3 Definition of jth level NDWT coefficients

For arbitrary sample size n, the NDWT wavelet and scaling coefficients are defined by

$\begin{array}{l}{\tilde{W}}_{j,t}={\displaystyle {\sum }_{k=0}^{K_j-1}{\tilde{h}}_{j,k}{X}_{t-kmodn}}\\ {\tilde{V}}_{j,t}={\displaystyle {\sum }_{k=0}^{K_j-1}{\tilde{g}}_{j,k}}{X}_{t-kmodn},t=\mathrm{0,1,...,}n-\mathrm{1,}\end{array}$ (2.7)

where $\left\{{\tilde{h}}_{j,k}:k=\mathrm{0,...,}{K}_j-1\right\}$ and $\left\{{\tilde{g}}_{j,k}:k=\mathrm{0,...,}{K}_j-1\right\}$ are, respectively, the NDWT wavelet and scaling filters defined via ${\tilde{h}}_{j,k}\equiv h_{j,k}/2^{j/2}$ and ${\tilde{g}}_{j,k}\equiv g_{j,k}/2^{j/2}$ from the wavelet {h _j,k } and scaling {g _j,k } filters of width K _j ≡ (2^j - 1)(K - 1) + 1.

The NDWT filters are modified at each scale by inserting zeros. That is, at each scale 2^{j - 1} zeros are entered between each K values of the NDWT filters $\left\{{\tilde{h}}_j\right\}$ and $\left\{{\tilde{g}}_j\right\}$, namely

$\begin{array}{l}{\tilde{h}}_0,\underset{2^{j-1}}{\underbrace{\mathrm{0,...,0}}},{\tilde{h}}_1,\underset{2^{j-1}}{\underbrace{\mathrm{0,...,0}}}\mathrm{,...,}{\tilde{h}}_{K-2},\underset{2^{j-1}}{\underbrace{\mathrm{0,...,0}}},{\tilde{h}}_{K-1}\\ {\tilde{g}}_0,\underset{2^{j-1}}{\underbrace{\mathrm{0,...,0}}},{\tilde{g}}_1,\underset{2^{j-1}}{\underbrace{\mathrm{0,...,0}}}\mathrm{,...,}{\tilde{g}}_{K-2},\underset{2^{j-1}}{\underbrace{\mathrm{0,...,0}}},{\tilde{g}}_{K-1},\end{array}$ (2.8)

that consist in apply an upsample of width 2^{j - 1} (K - 1) + 1. This process eliminate the downsample 5.5 ⁸.

2.3.1 The NDWT pyramidal algorithm

The NDWT wavelet coefficient ${\tilde{W}}_j$ and NDWT scaling coefficient ${\tilde{V}}_j$ of j level, presented in the equation (2.7), can still be calculated using an efficient algorithm, based on NDWT scaling coefficient ${\tilde{V}}_{j-1}$ of the j?1 level.The NDWT pyramidal algorithm is similar to the DWT one with the advantage that it does not produce a downsampling of wavelet and scaling coefficients and does not require the length of the signal to be a power of two. At each stage of the NDWT pyramidal algorithm the wavelet and scaling filters are upsampled, as in (2.8), so that when performing the convolving the signal with the filter, to obtain n coefficients at each scale of the algorithm.

The NDWT coefficients ${\tilde{W}}_j,{\tilde{V}}_j$ and ${\tilde{V}}_{j-1}$ are obtained by circularly filtering X _t with the respective periodized filters $\left\{{\tilde{h}}_{j,k}\right\},\left\{{\tilde{g}}_{j,k}\right\}$ e $\left\{{\tilde{g}}_{j-\mathrm{1,}k}\right\}$, as in equation (2.7).

Moreover, it can be shown (see section 5.5 ⁸) that it is possible to obtain ${\tilde{W}}_j$ and ${\tilde{V}}_j$ by filtering of ${\tilde{V}}_{j-1}$ by the following equation

$\begin{array}{l}{\tilde{W}}_{j,t}={\displaystyle {\sum }_{k=0}^{K-1}{\tilde{h}}_k{\tilde{V}}_{j-\mathrm{1,}t-2^{j-1}k\mathrm{mod}N}}\\ {\tilde{V}}_{j,t}={\displaystyle {\sum }_{k=0}^{K-1}{\tilde{g}}_k{\tilde{V}}_{j-\mathrm{1,}t-2^{j-1}k\mathrm{mod}N},t=\mathrm{0,1,...,}N-1}\end{array}$ (2.9)

These two equations in (2.9) constitute the NDWT pyramidal algorithm and also can be written using the N?N matrices ${\tilde{B}}_j$ and ${\tilde{A}}_j$, as

$\begin{array}{l}{\tilde{W}}_j={\tilde{B}}_j{\tilde{V}}_{j-\mathrm{1,}t}\\ {\tilde{V}}_j={\tilde{A}}_j{\tilde{V}}_{j-\mathrm{1,}t},\end{array}$

where the rows of ${\tilde{B}}_j$ consist of the wavelet filter $\left\{{\tilde{h}}_j\right\}$ upsampled width 2^{j - 1} (K - 1) + 1 and periodized to length n, with each row differing from its neighbors by circular shifts of one unit either forward ou backward. Likewise, the rows of ${\tilde{B}}_j$ are built based in the scaling filter $\left\{{\tilde{h}}_j\right\}$.

Note that defining ${\tilde{V}}_{\mathrm{0,}t}=X$, the equations in (2.9) produce wavelet coefficients ${\tilde{W}}_1$ and scaling coefficients ${\tilde{V}}_1$ of the first level of the NDWT.

The NDWT also allows to reconstruct ${\tilde{V}}_{j-1}$ from ${\tilde{W}}_j$ and ${\tilde{V}}_j$. The inverse NDWT can be calculated via inverse pyramidal algorithm described by the following equation

${\tilde{V}}_{j-\mathrm{1,}t}={\displaystyle \sum _{k=0}^{K-1}{\tilde{h}}_k{\tilde{W}}_{j,t+2^{j-1}k\mathrm{mod}N}+}{\displaystyle \sum _{k=0}^{K-1}{\tilde{g}}_k{\tilde{V}}_{j,t+2^{j-1}k\mathrm{mod}N},}t=\mathrm{0,1,...,}N-1.$

Or, in the matrix form

${\tilde{V}}_{j-1}={\tilde{A}}_j^t{\tilde{W}}_j+{\tilde{B}}_j^t{\tilde{V}}_j.$ (2.10)

Identifying ${\tilde{V}}_0\equiv X$ and applying (2.10) recursively until the stage J, can be expressed as

$X={\tilde{A}}_1^t{\tilde{W}}_1+{\tilde{B}}_1^t{\tilde{A}}_2^t{\tilde{W}}_2+{\tilde{B}}_1^t{\tilde{B}}_2^t{\tilde{A}}_3^t{\tilde{W}}_3+\mathrm{...}+{\tilde{B}}_1^t\mathrm{...}{\tilde{B}}_{J-1}^t{\tilde{A}}_J^t{\tilde{W}}_J+{\tilde{B}}_1^t\mathrm{...}{\tilde{B}}_{J-1}^t{\tilde{B}}_J^t{\tilde{V}}_J,$

considering the jth NDWT detail coefficient ${\tilde{D}}_j$ and J-th NDWT smooth coefficient ${\tilde{S}}_J$, namely

${\tilde{D}}_j={\tilde{B}}_1^t\mathrm{...}{\tilde{B}}_{j-1}^t{\tilde{A}}_j^t{\tilde{W}}_j\text{and}{\tilde{S}}_j={\tilde{B}}_1^t\mathrm{...}{\tilde{B}}_{J-1}^t{\tilde{B}}_J^t{\tilde{V}}_J$ (2.11)

and given a sample size n, can be expressed in NDWT additive decomposition by

$X={\displaystyle \sum _{j=1}^J{\tilde{D}}_j+{\tilde{S}}_J}$

where the energy of the NDWT decomposition can be obtained from

${\Vert X\Vert }^2={{\displaystyle \sum _{j=1}^J\Vert {\tilde{W}}_j\Vert }}^2+{\Vert {\tilde{V}}_J\Vert }^2$

for any integer J ³ 1.

The NDWT pyramidal algorithm is similar to the DWT one with the advantage that it does not produce a downsampling of wavelet and scaling coeffic? ients and does not require the length of the signal to be a power of two.

2.4 Signal analysis by wavelet periodogram

The wavelet periodogram is calculated from wavelet coefficients, namely

$I_{j,t}=\left|W_{j,t}\right|,$

decomposing the energy of the signal in multi-scales and with time-frequency location.

The Table (1) provides a scale-frequency interpretation of the wavelet periodogram.

Table 1:
Frequency and Period Intervals in each scale τ_j = 2^{j - 1}

Resolution Level	Scale	Frequency Interval	Period Interval
j	τj	$\left(\frac{1}{4{\tau }_j},\frac{1}{2{\tau }_j}\right]$	[2τj , 4τj )
1	τ1 = 1	$\left(\frac{1}{4},\frac{1}{2}\right]$	2 to 4
2	τ2 = 2	$\left(\frac{1}{8},\frac{1}{4}\right]$	4 to 8
3	τ3 = 4	$\left(\frac{1}{16},\frac{1}{8}\right]$	8 to 16
4	τ4 = 8	$\left(\frac{1}{32},\frac{1}{16}\right]$	16 to 32
5	τ5 = 16	$\left(\frac{1}{64},\frac{1}{32}\right]$	32 to 64
6	τ6 = 32	$\left(\frac{1}{128},\frac{1}{64}\right]$	64 to 128
7	τ7 = 64	$\left(\frac{1}{256},\frac{1}{128}\right]$	128 to 256
⋮	⋮	⋮	⋮

Each scale τ_j related to resolution level j corresponds to a period interval from 2τ_j to 4τ_j . The higher the resolution level j, the smoother the scale τ_j , representing the effects related to this scale have low frequency.

The periodogram provides a good description of the predominant frequencies contained in the signal and where the significant changes are located.

3 EXPERIMENTS

As already mentioned, the DWT can be seen as a filtering process followed by a downsampling, wherein only the odd or even elements of the sample are availed. Therefore, the result of the DWT depends on the choice of origin, i.e., the DWT is not shift-variant. To illustrate this fact, the DWT and NDWT were applied to both simulated and real data. In the case of simulated data, we generated signal s, illustrated in (b) of Figure (1), consisting of the sum of three periodic functions f, g and h of frequencies k = 12,21 and 5 (showed in (a) of Figure (1)), respectively, and a localized noise m added to the signal. In (c) of Figure (1), the same signal s is translated from 20 units forward.

Figure 1:
Simulated and shifted signals composed by periodic functions and a non-stationary behavior.

As all the components (f, g, h and m) of the simulated signal are known, we can evaluate the identification of these components after decomposition in periodogram analysis from both DWT and NDWT. Furthermore, after identifying specific behaviors depending on the frequencies or scales of interest in the periodogram, the reconstruction can be compared to the original signal, in order to assess the performance of these transforms in signal analysis. Furthermore, after identifying specific behaviors depending on the frequencies or scales of interest in the periodogram, the reconstruction can be compared to the original signal, allowing investigate which of the transforms (DWT or NDWT) generates a more perfect reconstruction of the signal. Thus, e.g. if we intend to identify and remove effects of frequency 12 from the signal s, after the decomposition, the scale that includes this frequency can be removed. Then, the coefficients of other scales can be reconstructed and compared to the original simulated signal (g + h + m) in the time domain.

To illustrate the shift-invariance of the NDWT also in the real data analysis, we analyzed the effects of ionospheric scintillation on GPS signals ². The ionospheric scintillation, caused by irregularities in the density of electrons present in the ionosphere, can weaken the signal received by the GPS receiver, causing degradation of positioning or even signal loss ³. The main indicative to investigate the ionospheric scintillation impact in GPS satellite signals quality is the S4 index. In (a) of Figure (2) we show the ionospheric scintillation S4 signal in a day of weak effect of scintillation. The time series of scintillation index presents daily gaps that occur when the satellite is invisible below the horizon, and daily periodic behavior, which has a shape of "U", when there are data. This scintillation S4 signal was shifted some units forward (in (b) of Figure (2)) and it also was analyzed in multi-scale by DWT and NDWT.

Figure 2:
Analyzed S4 index time series.

4 RESULTS AND DISCUSSIONS

4.1 Simulated data

The DWT and NDWT were applied to the simulated signal s and shifted one, in order to obtain a performance comparison of both transforms in signal analysis. The Figures (3) and (4) present the DWT and NDWT periodogram of the signals, respectively.

Figure 3:
DWT periodogram.

Figure 4:
NDWT periodogram.

The comparison of the obtained results in Figures (3) and (4) shows a better performance of the NDWT in the signal s analysis. Through DWT periodogram it is not easy to extract information of the signal. The NDWT periodogram, however, spells out the periodic functions on the scales of resolution levels from 2 to 4 in addition to identifying noise in the most refined scale, which is located in the signal between the instants 40 and 50. According to Table (1), the periodic functions f, g and h, respectively located on the scales of resolution levels 3,4 and 2, have frequencies in the bands from 8 to 16, 16 to 32 and 4 to 8, respectively, as expected.

In applications like signal analysis, which it is intended to investigate characteristics and patterns in the signal, it is necessary that the method be not sensitive to the origin, i.e., it be shift-invariant. We can observe in (b) of Figure (3) that the translation in data input results in a set of wavelet coefficients different of the data set obtained with the original data input (in (a) of Figure (3)), showing that the DWT is shift-variant. However, note in the Figure (4) that a translation of the input signal does not generate changes in the wavelet coefficients (in (b)) relative to the original signal (in (a)), unless than a translation. Functionality of the NDWT described herein motivates their use in the both stationary and non-stationary signal analysis.

Simulations were made in order to compare the performance of the DWT and NDWT in the reconstruction analysis. Each effect concerning to the frequencies of the functions f, g, h and m, which comprise the simulated signal s (in (b) of Figure (1)), were singly eliminated in the multi-scale decomposition. Then the signal was reconstructed and compared to the simulated signal without such frequency. The Figures (5) and (6) present the obtained results with the DWT and NDWT, respectively. It is important to be clear, that in (a) of these Figures, for example, the estimated signal f + h + m is the reconstruction from the complete signal s=f+g+h+m after decomposition, identification and removal of the scale of resolution level 4 (frequency band from 16 to 32) that contents the frequency of the signal g. Such estimated signal is compared with the original simulated signal f + h + m. This same procedure is performed in all other cases (subfigures (b), (c) and (d)).

Figure 5:
Reconstruction analysis by DWT. Simulated function and estimated function reconstructed from the complete signal s= f+g+h+m after decomposition, identification and removal of the scale that contains the frequency of interest.

Figure 6:
Reconstruction analysis by NDWT. Simulated function and estimated function reconstructed from the complete signal s= f+g+h+m after decomposition, identification and removal of the scale that contents the frequency of interest.

The Figures (5) and (6) show the reconstruction analysis by DWT and NDWT, respectively. Note by comparing these Figures that, in both instances, the estimation from the reconstructed signal was not perfect, as expected. This fact is justified by removing an entire bandwidth of the wavelet periodogram. Even so the estimated signal closely approximates the original signal, mainly from NDWT that presents better results than DWT.

In order to compare the quality of the reconstructions obtained from the DWT and NDWT (Figures (5) and (6)) we calculate the Mean Absolute Error (MAE) and Mean Squared Error (MSE) between the simulated signal and signal estimated from each wavelet decomposition. The results are presented in the Table 2.

Table 2:
Mean Absolute Error (MAE) and Mean Squared Error (MSE) of signal estimated from NDWT and DWT.

Estimated signal	MAE	MSE	MAE	MSE
f+h+m	0.59	0.67	0.99	1.50
g+h+m	0.83	1.38	1.31	2.75
f+g+m	0.43	0.49	0.99	1.43
f+h+g	0.52	1.43	0.73	1.76

Observe in the Table 2 that for all signal estimated from NDWT the MAE and MSE were lower than those estimated by DWT. It shows the NDWT superiority in relation to DWT in this analysis.

4.2 Real data

The DWT and NDWT were also applied in the S4 index time series ². The Figures (7) and (8) illustrate the DWT and NDWT wavelet periodogram of the S4 signal (in (a)) and shifted S4 signal (in (b)), respectively. Through them we compared the performance of both transforms in signal analysis.

Figure 7:
Multiscale analysis of S4 index time series by DWT.

Figure 8:
Multi-scale analysis of S4 index time series by NDWT.

In the S4 index time series multiscale analysis, the NDWT also present better results than DWT. The comparison of the periodograms produced by DWT (Figure (7)) and NDWT (Figure (8)), shows that differently of the set of wavelet coefficients generated by DWT, the NDWT coefficients of shifted signal (in (b) of figure) differ of the NDWT coefficients original S4 signal (in (a)) only by translation of coefficients, showing that NDWT is shift-invariant also in real data analysis.

5 CONCLUSIONS

Although NDWT has been underused in the literature, it can provide superior results in the signal analysis when compared to the DWT. The NDWT can also be calculated through Pyramidal Algorithm, similarly to DWT. The theory involved in the its calculation was discussed in this paper. The advantages of applying NDWT in signal analysis, such as, no sensitivity to origin choice (shift-invariant), detection of patterns and hidden information present in the time series and its time-frequency location, as well as perfect reconstruction were clearly showed from the simulated and real data. Hence, we encourage the use of NDWT instead of DWT when the aims is to analyze the coefficients or investigate effects in signals or time series.

ACKNOWLEDGMENTS

Thanks to Foundation for Research Support of the State of São Paulo - FAPESP and CAPES, for the financial support.

REFERENCES

• ¹
  M. Bašta. The MODWT analysis of the relationship between mortality and ambient temperature for Prague, Czech Republic, in Acta Oeconomica Pragensia, 19 (2011), 20-40.
• ²
  G. O. N. Brassarote, E. M. Souza & J. F. G. Monico. Multiscale Analysis of GPS Time Series from Non-decimated Wavelet to Investigate the Effects of Ionospheric Scintillation, in TEMA - Tend. Mat. Apl. Comput., 16(2) (2015), 119-130.
• ³
  R. S. Conker, M. B. El-arini, C. J. Hegarty & T. Hsiao. Modeling the Effects of Ionospheric Scintillation on GPS/Satellite-Based Augmentation System Availability, in Radio Science, 38 (2003), 1-23.
• ⁴
  M. González-Audícana, X. Otazu, O. Fors & A. Seco. Comparison between Mallat’s and the ’à trous’ discrete wavelet transform based algorithms for the fusion of multispectral and panchromatic images, in International Journal of Remote Sensing 26(3) (2005), 595-614.
• ⁵
  A. Kumar, L. K. Joshi, A. Pal & A. Shukla. MODWT based time scale decomposition analysis of BSE and NSE indexes financial time series, in International Journal of Mathematical Analysis, 5(27) (2011), 1343-1352.
• ⁶
  S. G. Mallat. A theory for multiresolution signal decomposition: the wavelet representation, in IEEE Trans. on Pattern Anal. and Mach.Intell., 11 (1989), 674-693.
• ⁷
  G. P. Nason & B. W. Silverman. The Stationary Wavelet Transform and Some Statistical Applications, in Wavelets and statistics - Springer Verlag, (1995), 281-281.
• ⁸
  D. B. Percival & A. T. Walden. Wavelet Methods for Time Series Analysis (Vol. 4), Cambridge university press, (2006).
• ⁹
  M. J. Shensa. Discrete Wavelet Transforms: The relationship of the ’à trous’ and Mallat algorithms, in Treizième Colloque Gretsi - Juan-Les-Pins, (1991).
• ¹⁰
  F. Buendia, A. M. Tarquis, G. Buendia & D. Andina. Feature Extraction Via Multiresolution MODWT Analysis in a Rainfall Forecast System, in Wmsci 2008: 12th World Multi-Conference on Systemics, Cybernetics and Informatics, (2008), 69-73.
• ¹¹
  C. Torrence & G. P. Compo. A Practical Guide to Wavelet Analysis in Bulletin of the American Meteorological Society, 79(1) (1998), 61-78.

Publication Dates

• Publication in this collection
  Jan-Apr 2018

History

• Received
  30 Jan 2017
• Accepted
  05 Mar 2018