Monday, June 11, 2018

matlab - Number of Daubechies coefficients



I am wondering about the correlation between input size and number of coefficients given by a discrete wavelet transform.


I am using Daubechies wavelets to describe a 1D function and I'm using PyWavelets to implement it (which is analogous to the MATLAB toolbox).


I started by implementing it using Haar wavelets, which gave correct results and I understand exactly how it works. Let's say my input function has 16 datapoints, if I use Haar, what I get from a multilevel decomposition (wavedec) is something like this (the number of shifts in brackets):


V1[1], W1[1], W2[2], W3[4], W4[8]

This is all well and good. The V1 gives me the scaling function and the W1-W5 wavelets of different scale and dilation. My problem is when I use the next Daubechies (referred to as 'db2' in the toolbox, which is called the D4), and I get


V1[6], W1[6], W2[9]

I lose all my intuition. I have no idea where 6, 6 and 9 come from, and they change depending on the level I specify (not even sure what it means to specify a level) and of course the input size. How can I predict the number of coefficients, and what are some good resources for gaining better understanding of why?


Thanks!



EDIT: Clarification on V and W:


$V_n$ usually denotes the span of a certain scaling function, $\phi$, i.e. $\{\phi_{n,k}\}$, where $k$ is the shift and $n$ the scaling. $W_n$ is the same except for the wavelet function, $\psi$. I might have abused the notation a bit by referring to the vectors of coefficients by V and W though.


EDIT2: Code


Here is the MATLAB code to produce the above:


>> [C, L] = wavedec(1:16, 4, 'db1'); L
L =
1 1 2 4 8 16
>> [C, L] = wavedec(1:16, 2, 'db2'); L
L =
6 6 9 16


I actually used PyWavelets, where it looked like this:


>>> import pywt
>>> map(len, pywt.wavedec(range(16), 'db1')) # defaults to level = 4
[1, 1, 2, 4, 8]
>>> map(len, pywt.wavedec(range(16), 'db2')) # defaults to level = 2
[6, 6, 9]

Answer



According to MATLAB documentation on wavedec,




The length of each filter is equal to $2N$. If n = length(s), the signals $F$ and $G$ are of length $n + 2N −1$ and the coefficients $cA_1$ and $cD_1$ are of length


$$\text{floor}\left( \frac{n-1}{2} \right) + N$$



Here, $n=16$ is the length of your signal, and $N=2$ is the Daubechies number.


Putting all that together, your detail coefficients at the second level should be of length


$$\text{floor}\left( \frac{16-1}{2} \right) + 2 = 7+2 = 9.$$


At the second level, you coefficients should be of length


$$\text{floor}\left( \frac{9-1}{2} \right) + 2 = 4+2 = 6.$$


If you're wondering why this must be the case, imagine the filtering-decimation procedure. Scaling and wavelet function for $\text{db}m$ wavelets are of length $2m$. When you convolve length $n$ signal with length $2m$ signal, you get back a signal of length $l_0=n+2m-1$. If you take every second sample of this resulting signal, you get back something of length $l_1 = \frac{n+2m-1}{2} = \frac{l_0}{2}$. Of course, if $l_0$ is odd ($n$ is even), then we cannot split the signal exactly into two parts. With clever mathematics that I will not go into here, you can show that this last non-paired coefficient is redundant anyway (carries no information that we don't yet know), so you can just omit it. Therefore, we will always have the resulting decimated signal of length


$$l_1 = \text{floor}\left( \frac{l_0}{2} \right) = \text{floor}\left( \frac{n+2N-1}{2} \right) = \text{floor}\left( \frac{n-1}{2} \right)+N.$$



No comments:

Post a Comment

periodic trends - Comparing radii in lithium, beryllium, magnesium, aluminium and sodium ions

Apparently the of last four, $\ce{Mg^2+}$ is closest in radius to $\ce{Li+}$. Is this true, and if so, why would a whole larger shell ($\ce{...