3
$\begingroup$

I am working on ECG signals, to eventually extract features in order to detect an arrhythmia and classify it. I am using Discrete Wavelet Transform with biorthogonal wavelet bior6.8 During my research, I came to know that wavelet transform is the convolution of the input processed signal with the daughter wavelets to get approximation and detail

$$X(a,b)=\frac{1}{\sqrt{a}}\int_{-\infty}^{\infty} \Psi\left( \frac{t-b}{a}\right)x(t)dt $$ where $a$ is scaling and $b$ is time.

I can't find the expression of the mother wavelet anywhere also when it came to practice DWT is usually presented as a filter bank of high pass and low pass filter The question is what is the difference between the different wavelets if it is always presented as a bank of filters filter bank

In my work, I used these two Butterworth high pass and low pass filters, but I still can't explain my choice, I read that Butterworth is the most used in signal processing and that it optimizes the frequency response in the passband, getting as much as you can from the wanted frequency Still, I have no arguments why I shouldn't use any others and not sure whether it is the correct way to implement bior 6.8 wavelet as I do not know any other way to implement wavelets and I would implement Daubechies or any other the same which do not make sense

from scipy.signal import filtfilt , butter

def butter_highpass(cutoff, fs, order=5):
    nyq = 0.5 * fs
    normal_cutoff = cutoff / nyq
    b, a = butter(order, normal_cutoff, btype='high', analog=False)
    return b, a

def butter_highpass_filter(data, cutoff, fs, order=5):
    b, a = butter_highpass(cutoff, fs, order=order)
    y = filtfilt(b, a, data,padlen=0)
    return y

from scipy.signal import butter, lfilter

def butter_lowpass(cutoff, fs, order=5):
    nyq = 0.5 * fs
    normal_cutoff = cutoff / nyq
    b, a = butter(order, normal_cutoff, btype='low', analog=False)
    return b, a

def butter_lowpass_filter(data, cutoff, fs, order=5):
    b, a = butter_lowpass(cutoff, fs, order=order)
    y = lfilter(b, a, data)
    return y

So basically and to resume, the question is what is the difference between wavelet transforms and how can we value that during implementation and then, the use of BF in this type of wavelet is it correct?

NB: the image inserted is from Wikipedia

$\endgroup$
  • 2
    $\begingroup$ So what is the question / problem? $\endgroup$ – Irreducible Mar 27 '18 at 13:14
  • $\begingroup$ I just edited the question , I want to know the difference between waveets in implementation $\endgroup$ – Abyr Mar 27 '18 at 13:29
0
$\begingroup$

The integral formula is the continuous wavelet transform, where you indeed convolve the data with a (mother) wavelet at all locations and scales. It is, by essence, redundant, and very mild conditions are required from continuous wavelet shapes to have both analysis and synthesis, and many shapes exist in analytic form.

If you want to discretize the above scheme in a critical way, ie to convert $N$ samples of a signal into $N$ coefficients, then you reach the realm of the discrete wavelet transform. Here, quite a limited choice of admissible wavelets exists (compared to the continuous case). And there, most continuous wavelets don't fit in the frame (pun intended). In other words, they cannot be discretized perfectly. The implementation technique was reversed: it was shown that the DWT could be implemented as the filter bank you mention. And from wavelet properties, filters properties where designed. Then, given a specific pair of admissible filters, wavelets are obtained through an iterative process, and generally, there is not closed form or formulae for such wavelets. If you are interested in more details, you can read Filters and Wavelets for Dyadic Analysis by D. S. G. Pollock:

Except in the polar cases of the Haar wavelet and the Shannon or sinc function wavelet, an explicit functional form for the wavelet or the scaling function is unlikely to be available. Nor is there, in most practical applications, an indispensable requirement to represent of these functions graphically

A whole DWT framework exists in Python: PyWavelets - Wavelet Transforms in Python. For instance, you can get filter coefficients for the biorthogonal 6.8 wavelet [EDIT: on the page wavelets.pybytes.com/wavelet/bior6.8 under Coefficients, click "Show values" (or "Copy")].

Conversely, if you take any pair of low/high pass filters, it will not necessarily become a genuine DWT filter bank, because they don't have the necessary "wavelet type" properties.

So, basically, to implement a DWT the standard way, you should pick filter coefficients tabulated in the literature, or listed at the PyWavelet property browser.

$\endgroup$
  • 1
    $\begingroup$ Can you please share with me a paper where I can find the filter coefficients for analysis and synthesis of filter bank corresponding to bior6.8 ? $\endgroup$ – Abyr Mar 27 '18 at 15:04
  • $\begingroup$ On the page wavelets.pybytes.com/wavelet/bior6.8 under Coefficients, click ""Show values" (or "Copy") $\endgroup$ – Laurent Duval Mar 27 '18 at 15:11
3
$\begingroup$

You've mentioned Butterworth filters for doing the wavelet analysis using bior6.8. If you want to perform the Discrete Wavelet Transform using some specific wavelet, then you must use its Perfect Reconstruction Filter Bank. Each wavelet function has its associated set of filter values for decomposition and reconstruction - they are calculated from the Mother and Daughter wavelet.

Since you were using Python, I am going to give you an example using PyWavelets. It is as simple as creating a wavelet object and asking for the filter bank. Here is a plot of the filter values:

enter image description here

Since these are FIR filter coefficients, then we can calculate the amplitude response by taking the DFT of the coefficient values:

enter image description here


Values of filter coefficients:

dec_lo = [0.0,
 0.0019088317364812906,
 -0.0019142861290887667,
 -0.016990639867602342,
 0.01193456527972926,
 0.04973290349094079,
 -0.07726317316720414,
 -0.09405920349573646,
 0.4207962846098268,
 0.8259229974584023,
 0.4207962846098268,
 -0.09405920349573646,
 -0.07726317316720414,
 0.04973290349094079,
 0.01193456527972926,
 -0.016990639867602342,
 -0.0019142861290887667,
 0.0019088317364812906]

dec_hi = [-0.0,
 0.0,
 -0.0,
 0.014426282505624435,
 -0.014467504896790148,
 -0.07872200106262882,
 0.04036797903033992,
 0.41784910915027457,
 -0.7589077294536541,
 0.41784910915027457,
 0.04036797903033992,
 -0.07872200106262882,
 -0.014467504896790148,
 0.014426282505624435,
 -0.0,
 0.0,
 -0.0,
 0.0]

rec_lo = [0.0,
 0.0,
 0.0,
 0.014426282505624435,
 0.014467504896790148,
 -0.07872200106262882,
 -0.04036797903033992,
 0.41784910915027457,
 0.7589077294536541,
 0.41784910915027457,
 -0.04036797903033992,
 -0.07872200106262882,
 0.014467504896790148,
 0.014426282505624435,
 0.0,
 0.0,
 0.0,
 0.0]

rec_hi = [0.0,
 -0.0019088317364812906,
 -0.0019142861290887667,
 0.016990639867602342,
 0.01193456527972926,
 -0.04973290349094079,
 -0.07726317316720414,
 0.09405920349573646,
 0.4207962846098268,
 -0.8259229974584023,
 0.4207962846098268,
 0.09405920349573646,
 -0.07726317316720414,
 -0.04973290349094079,
 0.01193456527972926,
 0.016990639867602342,
 -0.0019142861290887667,
 -0.0019088317364812906]

A complementary code which might be useful for you to investigate other wavelets:

import pywt
import numpy as np
import matplotlib.pyplot as plt

# Create wavelet and extract the filters
wavelet_name = 'bior6.8'
wavelet = pywt.Wavelet(wavelet_name)
dec_lo, dec_hi, rec_lo, rec_hi = wavelet.filter_bank

# Filter coefficients
plt.figure()
plt.subplot(221)
plt.stem(dec_lo)
plt.grid()
plt.title('{} low-pass decomposition filter'.format(wavelet_name))
plt.subplot(222)
plt.stem(dec_hi)
plt.grid()
plt.title('{} high-pass decomposition filter'.format(wavelet_name))
plt.subplot(223)
plt.stem(rec_lo)
plt.grid()
plt.title('{} low-pass reconstruction filter'.format(wavelet_name))
plt.subplot(224)
plt.stem(rec_hi)
plt.grid()
plt.title('{} high-pass reconstruction filter'.format(wavelet_name))

# Frequency responses
dec_lo_fr = np.abs(np.fft.rfft(dec_lo, 128))
dec_hi_fr = np.abs(np.fft.rfft(dec_hi, 128))
rec_lo_fr = np.abs(np.fft.rfft(rec_lo, 128))
rec_hi_fr = np.abs(np.fft.rfft(rec_hi, 128))

plt.figure()
plt.subplot(211)
plt.plot(dec_lo_fr, label='Low-pass')
plt.hold(True)
plt.plot(dec_hi_fr, label='High-pass')
plt.grid()
plt.legend()
plt.title('Frequency responses of {} decomposition filters'.format(wavelet_name))
plt.subplot(212)
plt.plot(rec_lo_fr, label='Low-pass')
plt.hold(True)
plt.plot(rec_hi_fr, label='High-pass')
plt.grid()
plt.legend()
plt.title('Frequency responses of {} reconstruction filters'.format(wavelet_name))

plt.show()

Regarding the actual filtering process you can either use scipy.signal.lfilter or np.convolve (x being the input signal):

y = lfilter(dec_lo, 1, x)

or

y = convolve(x, dec_lo)

I would suggest using PyWavelets to do all of that for you, instead of re-inventing the wheel. However, If you really want to implement a single level of DWT on your own then I second the numpy option. Here is a complete example:

t = np.linspace(0, 1.0, 128)
x = np.sin(2*np.pi*10*t)

# Perform manual DWT and decimate
cA = np.convolve(x, dec_lo)[1::2]
cD = np.convolve(x, dec_hi)[1::2]

plt.figure()
plt.subplot(211)
plt.plot(cA)
plt.grid()
plt.title('Approximation coefficients')

plt.subplot(212)
plt.plot(cD)
plt.grid()
plt.title('Detail coefficients')

plt.show()

Which gives:

enter image description here

$\endgroup$
  • $\begingroup$ So knowing the filters coefficients can I implement them ? $\endgroup$ – Abyr Mar 27 '18 at 14:47
  • $\begingroup$ Do you need anything else? Everything is there. Is anything unclear to you? $\endgroup$ – jojek Mar 27 '18 at 15:07
  • $\begingroup$ Well I'm new to signal processing , so I'm not sure how to implement a filter knowing the coefficients $\endgroup$ – Abyr Mar 27 '18 at 15:21
  • $\begingroup$ Please see my answer updated. You've got everything that you need to perform the analysis. Just don't use filtfilt. $\endgroup$ – jojek Mar 27 '18 at 15:26
  • 1
    $\begingroup$ @Abyr ... your question was of theoretical nature and I would say it was answered more than clear. Your follow up questions are neither related to wavelets nor signal processing. They are targeting programming, here I highly recommend you to read the help of these functions and google some tutorials. This way you will learn by far more ... $\endgroup$ – Irreducible Mar 28 '18 at 5:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.