A FAST WAVELET TRANSFORM-DOMAIN LMS ALGORITHM

S. Attallah and M. Najim

Equipe Signal et Image and GDR 134-CNRS
ENSERB. AV. du Dr. Albert Schweitzer
BP 99 - 33402 Talence Cedex
FRANCE
e-mail: [email protected]
ABSTRACT Second, we generalize the method to include longer
wavelet filters and bigger decimation factors. For both
This paper presents a Fast Wavelet Transform-Domain cases, we will give explicit formulae for the computation
LMS (FWTDLMS) algorithm whose computational of the complexity of the wavelet transform. Finally, we
complexity is very low as compared to what has been turn our attention to the self-orthogonalizing algorithm
reported so far, because it is completely independent from where, taking benifit from redundancy as mentioned above,
the adaptive filter order. In our method, the reduction in we will show how to estimate the energy of the wavelet
the computational load has been achieved by exploiting coefficients at a very low computational cost.
the redundancy in the computation of the wavelet
coefficients as the data window slides on by one sample at 2. TWO-BAND MAXIMALLY DECIMATED
each new iteration. As a result, the complexity of the FILTER BANK
wavelet transform is reduced from a") where m is the
wavelet filter length and N is the adaptive filter order to In this section, we consider the case of a 2-band filter bank
O(m(m+M)/2) where M (M < m) is the number of (fig. 1, M=2). Figure 2 shows how the data window slides
subbands in the maximally decimated filter bank. In on at each new iteration. For simplicity, only one subband
addition,we have introduced a recursion in the estimation and a short wavelet filter are considered. In the following,
of the energy of the wavelet coefficients. This complexity we consider that N, m and M are positive integers such
becomes then 0(2N+m(m+M)D+m+4M), which is that M < m < N and N and m are multiples of M. In
comparable to time-domain LMS complexity when short addition, the orthogonal wavelet transform matrix is
wavelet filters (m<<N) are used. constructed from the wavelet filters following [2], i.e., the
convolution boundary effects are handled by periodizing
1. INTRODUCTION the data vector so as to have an orthogonal matrix within
the interval of time whose length is equal to order N. In
The flexibility of the orthogonal wavelet transform, which figure 2, we show how data (x(n), x(n-l), ...) and the
lies in the possibility of designing filters very suitable to wavelet coefficients (ch(O), ch(l), ...) are arranged:
the application at hand, inspired researchers to use it as an periodization is carried out only at one side (see fig. 2) by
alternative to DCT and FFT for improving the wrapping around the wavelet filter coefficients at the
convergence of LMS algorithm when dealing with colored discrete times n, n+l, and so on. We can easily notice that
input signals 111, [2]. However, the computation of the even numbered blocks, i.e., block n and block (n+2)
wavelet transform introduces a complexity of order (mN), (respectively odd numbered blocks, i.e., block(n+l) and
where m and N are, respectively, the wavelet filter length block(n+3)) have at least (N-m)/2 wavelet coefficients
and the adaptive filter order. It is clear that for filters even which are the same. If these wavelet coefficients are stored,
as small as 4,8 and 16 coefficients the complexity of the then at each iteration, we would need to compute only m/2
wavelet transform by itself is 4N, 8N and 16N, wavelet coefficients. The computational complexity for
respectively. This is generally not acceptable. In this the two subbands is therefore
paper, we present a new method for reducing considerably
(figure 2) this complexity. Our method takes into
consideration the redundancy in the calculation of the
wavelet coefficients as the data window slides on by one Figure 2 also shows that,besides the number of wavelet
sample at each iteration. Moreover, we will show that the coefficients common to even numbered blocks
estimation of the energy of the wavelet coefficients can be (respectively odd numbered blocks), there is a certain
computed recursively which reduces still further the whole number of intermediate results that can be stored and used
complexity of the wavelet transform-domain LMS to advantage as well. thereby the computational
algorithm. The rest of the paper is organized as follows: complexity for the two subbands becomes (m being a
first we consider the case of a 2-band filter bank (M=2). multiple of 2)

0-7803-3192-3/96 $5.0001996 IEEE 1343

 C2 = 2 [m+(m-2)+(m-4)+...+(m+(m-2))] (2)
It is very easy to show that this sum has a mathematical The estimated signal is then given by
compact form (see appendix) which is given by

~2 = m (m + 2) (3)
2 The error to be minimized is defined as:
Therefore for a Two-band tree, wether it is dyadic or
uniform, the computational complexity is O(m(m+2)/2). e(@= d(n) - ycst(n) (11)

3. M-BAND MAXIMALLY DECIMATED The adaptive filter coefficients are updated using the
FILTER BANK following equations [31

Results obtained above can be readily extended to M-band

maximally decimated filter banks.For each subband and at
each iteration, each pair of block n and block(n+M) will
have in common (N-M)/M wavelet coefficients. So it will
remain only Mm/M wavelet coefficients to be computed
for the whole M-band filter bank and which can be reduced where &o(n) and &(n) are the energy estimates of the
still further by storing some intermediate results (at the wavelet coefficients at the nth iteration. These are obtained
previous iteration) such that iteratively using the following two equations

CM= M [m+(m-M)+(m-2M)+...+( m-(m-M))] (4) &An)= p &(n-1) + (1 - p) pho(n>12 (14)

After development and rearrangement (m being a multiple &(n) = p &(n-1) + (1 - p) Dhl<n>12 (15)
of M), we get (see appendix) and
CM=m(m + M) (3 Bho(n>12= do(n) + bo(n-1) +...+ &o(n-$kl) (16)
2
We can notice from figure 3 that this computational a and p are positive constants less than 1.
complexity increases very slowly as the decimation factor We know from the previous discussion that in eqt (16) (N-
or the number of subbands increases. Therefore, our m)D wavelet coefficients have been already used to
method (O(m(m+M)/2) outperforms by far the direct compute kho(n-2j2 . So, if we have stored this partial
method @M) which is O(mN) [l]. result, then to compute pho(n12, we would need to
perform only m/2 multiplications per sample for each
4. WAVELET TRANSFORM-DOMAIN LMS subband. Thereby our method presents a computational
ALGORITHM complexity for the whole algorithm (wavelet transform +
4-1. FWTDLMS BASED ON A 2 - B A N D coefficients updating + estimation of the energy of the
MAXIMALLY DECIMATED FILTER BANK wavelet coeficients...) which is equal to

Let Z(n) be the wavelet transformed data column vector of C-ms = 2N + m ( m + 2) + m +8 (17)
length N, that is 2
E.g.: m=4, CFWTDLMS = 2N + 24 (comparable to the
Z(n) = TwX(n) (6) normalized LMS complexity). Whereas the direct method
Where @M) presents the following computational complexity
X(n) = [x(n) x(n-1) ....x(n-N+l)]' (7)
and (fig. 1)
z(n) = Pho(n) Zfil(n)f (8) = 7N + 8 (m = 4 is considered to
E.g.: m=4, C~irea-~eth~d
be very short !)
Data are supposed to be real, the superscript t to denote
vector transposition and Tw to represent the orthogonal For relatively long adaptive FIR filters, our method has a
wavelet transform matrix. Now, let g(n) be the adaptive computational complexity very comparable to that of a
filter column vector that we split down into two subfilters simple time-domain normalized LMS. Herebelow, we give
of length N/2 each such that a table (fig. 4) where we show the improvement of our

1344
method as compared to the direct method for different k
wavelet filter lengths and different values of N. We define C = C (m - i M)
i=n
the reduction in the complexity as compared to the direct
method0M) as
= (k-+ 1) m - M (1 + 2 + 3 + ...+ k) (24)
= (k+ 1 ) m - M- k ( k + 1 ) (25)
n
L
Red= 1 - CFWIDLMS
~Dirra-Method
= (k+ 1) (m - kM)
2
where CFWTDLMS is the complexity of our fast algorithm NOW replacing k in eqt (26) by m/M drawn from eqt (21)
and C~irca-~ahod
is the complexity of the same algorithm 1- to
based on the direct method for computing the wavelet C = m (m + M) (27)
transform. 2M
4-2 FWI"I'LMS BASED ON A M-BAND For a filter bank of M subbands, eqt (27) leads to eqt (5),
MAXIMALLY DECIMATED FILTER BANK i.e.,
CM= (m + M)
Extension of the above algorithm to M-bands is 2
straightforward. Following thesame way as before, it is Now if we replace M by 2 we get the result given by
very easy to show that for a M-band maximally decimated eqt (3).
filter bank, the computational complexity of the
FWTDLMS is given by: 7. REFERENCES

cmLm= 2~ + LE (m + M) + m + 4~ (20) [l] S. Hosur and A. H. Tewfik, "Wavelet Transform

2 Domain LMS Algorithm," Proc. ICASSP-93, pp. 508-
510, April 1993, Minneapolis, Minnesota, USA.
against
[2] S. Attallah and M. Najim, "On The Convergence
= 3N + mN + 4M
C~i~t-~cthod (21) Enhancement of The Wavelet Transform Based LMS,"
Proc. ICASSP-95, pp. 973-976 May 1995, Detroit,
for the direct-method based algorithm. Michigan, USA.
5. CONCLUSION [3] K. F. Wan and P.C. Ching, "On the Optimality of
Convergence Behaviour for Transform-Domain Split-Path
In this paper, we have presented a new fast wavelet Adaptive Filter" Proc. ICASSP-94, pp. 421424, April
transform-domain algorithm which takes into account the 1994, Adelaide, South Australia.
redundancy in the calculation of the wavelet transform, as
the data window slides on by one sample at each iteration, [41 J. C. Lee and C. K. Un, "Performance of Transform-
to reduce the computational load considerably. For large Domain LMS Adaptive Digital Filters," IEEE Trans.
values of N @ greater
i or equal to 64),the computational Acoust., Speech, Signal Processing, vol. ASSP-34, pp.
load is r e d u d by about 70%. For N=512, this reduction 499-510, June 1986.
can attain 90% (e.g. m=24, M=2 or 4). The reduction in
complexity does not affect the convergence performance of
the FWTDLMS which is similar to the one we presented

VZh
in [2]. The proposed method is computationally more
efficient than other transform based algorithms such as
DCT. because the latter has a computational complexity
0((6N+l)+NTR) where NTR = N log2(N) - (3/2)N +4,
[3], which is clearly higher than our method's one.
6. APPENDIX I I
In this appendix, we proove the results given by eqt (3)
and eqt (5). Since m is by hypothesis a multiple of M,
then we should have
I 1 I

zhM-l
m=kM (22)
- 1:used
Fin. M-band maximally decimated filter bank that is
to generate orthognal M-band wavelets
where k is a positive integer. Eqt (4) can be rewritten as

1345
 nth iteration

x(nt1) x(n) "1) x(n-2) x(n-3) x(n-4) x(n-5) x(n4)

Ch(1) Ch(2)
E$] Ch(2) m(3)
Ch(0) Ch(1) Ch(2) Ch(3)
(n+l)lh iteration QIM CYl) Q I M Ch(3)

(to be computed) I 1 Wawbtblockcoeficients

n&
mmmm to
block(n&\
0

Fig. 2: Redundancy in the calculation of the wavelet coefficients as the data window
slides on by one sample at each iteration (only one subband is shown).

Filter length m Fig. 4: FWTDLMS computationalcomplexity reduction

(in %) as compared to the direct method based algorithm
Fig. 3: Wavelet transform Computational complexity -
(Red = 1 C ~ I h f S)
comparaison between the proposed method and the direct CDirea-Mahod
method (DM)

1346