離散小波變換

Template:NoteTA 在数值分析和泛函分析领域中，离散小波变换（Discrete Wavelet Transform，DWT）是小波被离散采样的小波变换。与其他小波变换一样，它与傅里叶变换相比的一个关键优势是时间分辨率：它既能捕获频率信息，又能捕获位置（时间上的位置）信息。

第一個離散小波變換由匈牙利數學家哈尔發明，離散小波轉換顧名思義就是離散的輸入以及離散的輸出，但是這裡並沒有一個簡單而明確的公式來表示輸入及輸出的關係，只能以階層式架構來表示。

定義

• 首先我們定義一些需要用到的信號及濾波器。
• x[n]：離散的輸入信號，長度為N。
• $g[n]$：低通濾波器（low pass filter），可以將輸入信號的高頻部份濾掉而輸出低頻部份。
• $h[n]$：高通濾波器（high pass filter），與低通濾波器相反，濾掉低頻部份而輸出高頻部份。
• $\downarrow$ Q：降采样濾波器（downsampling filter），如果以x[n]作為输入，則輸出y[n]=x[Qn]。此處舉例Q=2。

舉例說明:

清楚規定以上符號之後，便可以利用階層架構來介紹如何將一個離散信號作離散小波轉換：

架構中的第1層(1st stage)

$x_{1,L}[n]=\sum _{k=0}^{K-1}x[2n-k]g[k]$

$x_{1,H}[n]=\sum _{k=0}^{K-1}x[2n-k]h[k]$

架構中的第2層(2nd stage)

$x_{2,L}[n]=\sum _{k=0}^{K-1}{x}_{1,L}[2n-k]g[k]$

$x_{2,H}[n]=\sum _{k=0}^{K-1}{x}_{1,L}[2n-k]h[k]$

可繼續延伸

$⋮$

$⋮$

架構中的第$\alpha$層($\alpha -th$ stage)

$x_{\alpha ,L}[n]=\sum _{k=0}^{K-1}{x}_{\alpha -1,L}[2n-k]g[k]$

$x_{\alpha ,H}[n]=\sum _{k=0}^{K-1}{x}_{\alpha -1,L}[2n-k]h[k]$

注意：若輸入信號$x[n]$的長度是N，則第$\alpha$層中的$x_{\alpha ,L}[n]$及$x_{\alpha ,H}[n]$的長度為$\frac{N}{2^{\alpha }}$

二維離散小波轉換

此時的輸入信號變成$x[m,n]$，而轉換過程變得更複雜，說明如下：

首先對n方向作高通、低通以及降頻的處理

$v_{1,L}[m,n]=\sum _{k=0}^{K-1}x[m,2n-k]g[k]$

$v_{1,H}[m,n]=\sum _{k=0}^{K-1}x[m,2n-k]h[k]$

接著對$v_{1,L}[m,n]$與$v_{1,H}[m,n]$延著m方向作高低通及降頻動作

$x_{1,LL}[m,n]=\sum _{k=0}^{K-1}{v}_{1,L}[2m-k,n]g[k]$

$x_{1,HL}[m,n]=\sum _{k=0}^{K-1}{v}_{1,L}[2m-k,n]h[k]$

$x_{1,LH}[m,n]=\sum _{k=0}^{K-1}{v}_{1,H}[2m-k,n]g[k]$

$x_{1,HH}[m,n]=\sum _{k=0}^{K-1}{v}_{1,H}[2m-k,n]h[k]$

經過(1)(2)兩個步驟才算完成2-D DWT的一個stage。

實際範例

以下根據上述2-D DWT的步驟，對一張影像作二維離散小波轉換(2D Discrete Wavelet Transform)

原始影像

2D DWT的結果

$\Rightarrow$

複雜度(Complexity)

在討論複雜度之前，先做一些定義，當x[n]*y[n]時，x[n]之長度為N，y[n]之長度為L： $\Rightarrow IDF{T}_{N+L-1}[DF{T}_{N+L-1}(x[n])DF{T}_{N+L-1}(y[n])]$

其中，

$IDF{T}_{N+L-1}$為(N+L-1)點離散傅立葉反轉換(inverse discrete Fourier transform)

$DF{T}_{N+L-1}$為(N+L-1)點離散傅立葉轉換(discrete Fourier transform)

(1)一維離散小波轉換之複雜度(沒有分段卷积(sectioned convolution))：

$\frac{3}{2}(N+L-1){\log }_2(N+L-1)\approx \frac{3}{2}N{\log }_{2}N$

(2)當 N >>> L 時，使用 “分段卷积(sectioned convolution)”的技巧：

將x[n]切成很多段，每段長度為${\mathbb{N}}_1$，總共會有$S=\frac{N}{N_1}$段，其中$N>N_1>>L$。

則

$x[n]*g[n]=x_1[n]*g[n]+x_2[n]*g[n]+...+x_s[n]*g[n]$

$x[n]*h[n]=x_1[n]*h[n]+x_2[n]*h[n]+...+x_s[n]*h[n]$

複雜度為：

$\begin{array}{cc}\frac{3}{2}S(N_1+L-1){\log }_2(N_1+L-1)& \approx \frac{3}{2}S{N}_1{\log }_2(N_1+L-1)\\ & \approx \frac{3}{2}N{\log }_2(N_1+L-1)\\ & \approx \frac{3}{2}N{\log }_2{N}_1\\ \end{array}$

在這裡要注意的是，當N>>L時，一維離散小波轉換之複雜度是呈線性的(隨N)，$𝑂\mathit{(}𝑁\mathit{)}$。

(3)多層(Multiple stages )的情況下：

1.若$x_{a,H}[n]$不再分解時：

$\begin{array}{cc}Complexity& \approx (N+\frac{N}{2}+\frac{N}{4}+\frac{N}{8}+...+2)\frac{3}{2}{\log }_2{N}_1\\ & =(2N-2)\frac{3}{2}{\log }_2{N}_1\\ & \approx 3N{\log }_2{N}_1\\ \end{array}$

2.若$x_{a,H}[n]$也細分時：

$\begin{array}{cc}Complexity& \approx (N+2\frac{N}{2}+4\frac{N}{4}+8\frac{N}{8}+...+\frac{N}{2}2)\frac{3}{2}{\log }_2{N}_1\\ & =(N{\log }_{2}N)\frac{3}{2}{\log }_2{N}_1\\ \end{array}$

(4)二維離散小波轉換之複雜度(沒有分段卷积(sectioned convolution))： $\Rightarrow M\frac{3}{2}(N+L-1){\log }_2(N+L-1)+(N+L-1)\frac{3}{2}(M+L-1){\log }_2(M+L-1)$

上式中，第一部分需要M個一維離散小波轉換並且每個一維離散小波轉換的輸入有N個點；第二部分需要N+L-1個一維離散小波轉換並且每個一維離散小波轉換的輸入有M個點。

$\begin{array}{cc}Complexity& \approx \frac{3}{2}MN{\log }_{2}N+\frac{3}{2}MN{\log }_{2}M\\ & =\frac{3}{2}MN({\log }_{2}N+{\log }_{2}M)\\ & =\frac{3}{2}MN{\log }_{2}MN\\ \end{array}$

(5)二維離散小波轉換之複雜度，使用 “分段卷积(sectioned convolution)”的技巧：

假設原始尺寸為$M\times N$，則每一小部分的尺寸為$M_1\times N_1$

$\begin{array}{cc}Complexity& \approx \frac{MN}{M_1{N}_1}\frac{3}{2}{M}_1{N}_1{\log }_2{M}_1{N}_1\\ & =\frac{3}{2}MN{\log }_2{M}_1{N}_1\\ \end{array}$

所以若是使用分段摺積，則二維離散小波轉換之複雜度是呈線性的(隨MN)，$𝑂\mathit{(}𝑀𝑁\mathit{)}$。

(6)多層(Multiple stages )與二維的情況下：

首先x[m,n]的尺寸為$M\times N$，

1.若$x_{a,H_1}[n],x_{a,H_2}[n],x_{a,H_3}[n]$不細分，只細分$x_{a,L}[n]$時，總複雜度為：

$\begin{array}{cc}totalcomplexity& =(MN+\frac{MN}{4}+\frac{MN}{16}+...)\frac{3}{2}{\log }_2{M}_1{N}_1\\ & \approx \frac{4}{3}MN\frac{3}{2}{\log }_2{M}_1{N}_1\\ & =2MN{\log }_2{M}_1{N}_1\\ \end{array}$

2.若$x_{a,H_1}[n],x_{a,H_2}[n],x_{a,H_3}[n]$也細分時，總複雜度為：

$\begin{array}{cc}totalcomplexity& =(MN+4\frac{M}{2}\frac{N}{2}+16\frac{M}{4}\frac{N}{4}+...)\frac{3}{2}{\log }_2{M}_1{N}_1\\ & =[MN{\log }_2(min(M,N))]\frac{3}{2}{\log }_2{M}_1{N}_1\\ \end{array}$

重建(Reconstruction)

使用離散小波轉換，將訊號個別通過一個低通濾波器和一個高通濾波器，得到訊號的高低頻成分，而在重建(Template:Link-en)原始訊號的過程，也就是離散小波的逆轉換(Inverse Discrete Wavelet Transform. IDWT)，直觀而言，我們僅是需要將離散小波轉換做重建濾波即可得到原始輸入信號，以下將推導重建濾波器，也就是IDWT高低通濾波器的構成要件，以及如何來重建原始信號。

重建過程如下：

使用Z轉換：

$X(z)=\sum x(n)z^{-n}$

• DWT低通濾波器 $g(n)$ 的Z轉換為 $G(z)$ ，DWT高通濾波器$h(n)$的Z轉換為$H(z)$

• 信號$x(n)$通過濾波器 $g(n)$ 後，Z轉換為 $X(z)$$G(z)$ ，信號$x(n)$通過濾波器 $h(n)$ 後，Z轉換為 $X(z)$$H(z)$

• 降低採樣數(downsample)2倍後，

$X_{1,L}(z)=\frac{1}{2}[X(z^{\frac{1}{2}})G(z^{\frac{1}{2}})+X(-z^{\frac{1}{2}})G(-z^{\frac{1}{2}})]$

$X_{1,H}(z)=\frac{1}{2}[X(z^{\frac{1}{2}})H(z^{\frac{1}{2}})+X(-z^{\frac{1}{2}})H(-z^{\frac{1}{2}})]$

• 升頻(interpolation)2倍後，再通過IDWT的低通重建濾波器 $g_1(n)$ ，

$X_0(z)=\frac{1}{2}[X(z)G(z)+X(-z)G(-z)]G_1(z)+\frac{1}{2}[X(z)H(z)+X(-z)H(-z)]H_1(z)$

$=\frac{1}{2}[G(z)G_1(z)+H(z)H_1(z)]X(z)+\frac{1}{2}[G(-z)G_1(z)+H(-z)H_1(z)]X(-z)$

若要完整重建，則$X_0(z)=X(z)$

條件1：$G(z)G_1(z)+H(z)H_1(z)=2$

條件2：$G(-z)G_1(z)+H(-z)H_1(z)=0$

因此，在設計高低通重建濾波器時，需要考慮上述條件，寫成矩陣形式如下：

$(\genfrac{}{}{0ex}{}{G_1(z)}{H_1(z)})=\frac{2}{det(H_m(z))}(\genfrac{}{}{0ex}{}{H(-z)}{-G(-z)})$

其中 $det(H_m(z))=G(z)H(-z)-H(z)G(-z)$

離散小波轉換(DWT)設計

四大條件：

1.DWT通濾波器 $g(n)$，$h(n)$必須要是有限長度。

2.滿足$h(n)$是高通濾波器(high pass filter)，$g(n)$是低通濾波器(low pass filter)。

3.滿足完整重建要條件，$(\genfrac{}{}{0ex}{}{G_1(z)}{H_1(z)})=\frac{2}{det(H_m(z))}(\genfrac{}{}{0ex}{}{H(-z)}{-G(-z)})$，其中 $det(H_m(z))=G(z)H(-z)-H(z)G(-z)$

4.若$g(n)$，$h(n)$為有限長度，則$det(H_m(z))=G(z)H(-z)-H(z)G(-z)=\alpha z^k$ ，且 $k$ 為奇數。

>第4點較難達成，是DWT設計的核心

*為什麼k是奇數？

假設k为偶数，

z=-1

$det(H_m(-1))=G(-1)H(1)-H(-1)G(1)=\alpha (-1{)}^k=1$

z=1

$det(H_m(1))=G(1)H(-1)-H(1)G(-1)=\alpha (1{)}^k=1$

$G(1)H(-1)-H(1)G(-1)=G(-1)H(1)-H(-1)G(1)$

$H(1)G(-1)=G(1)H(-1)$

代回

$det(H_m(-1))=G(-1)H(1)-H(-1)G(1)=0$

顯然出現矛盾。

所以k必須為奇数。

以下介紹兩種完美重建的DWT濾波器：

1.正交镜象滤波器（Quadrature Mirror Filter，QMF）

• $H(z)=G(-z)$

$⟹$ $h(n)=(-1{)}^{n}g(n)$

• $G_1(z)=G(z)z^{-k}$

$⟹$ $g_1(n)=g(n-k)$

• $H_1(z)=-G(-z)z^{-k}$

$⟹$ $h_1(n)=(-1{)}^{n-k+1}g(n-k)$

$det(H_m(z))=G(z)H(-z)-H(z)G(-z)=2z^k$，且 $k$ 為奇數。

2.单位正交小波（Orthonormal Wavelet）

• $G(z)=G_1(z^{-1})$

$⟹$ $g(n)=g_1(-n)$

• $H(z)=-z^k{G}_1(-z)$

$⟹$ $h(n)=(-1{)}^n{g}_1(n+k)$

• $H_1(z)=-z^k{G}_1(-z^{-1})$

$⟹$ $h_1(n)=(-1{)}^n{g}_1(-n+k)$

$det(H_m(z))=G(z)H(-z)-H(z)G(-z)=2z^k$，且 $k$ 為奇數。

多數小波屬於单位正交小波。

其他應用

壓縮、去除雜訊：使用低通濾波器，將小波轉換的高頻濾掉，即保留$x_{1,LL}[m,n]$而將其他部分捨棄。

• 邊緣偵測：使用高通濾波器，將小波的低頻濾掉，即保留$x_{1,HL}[m,n]$或$x_{1,LH}[m,n]$而捨棄其他部分。
• R语言小波分析wavelet Template:Wayback
• 作為 JPEG2000 的內部架構
• 模式辨認：由於可以利用低頻的部分得到原圖的縮略版，加上模式通常為整體的特性，藉由在縮略圖上進行工作，小波轉換可以有效減少尋找模式與比對模式的運算時間
• 濾波器設計：小波轉換保留部分時間資訊，可以據此資訊加上訊號的強度資訊，保留特定時點的資訊而同時去除雜訊

同時參閱

• 小波分析
• 連續小波轉換
• 哈爾小波轉換
• 信號處理

參考

• Jian-Jiun Ding, Time frequency analysis and wavelet transform class note,the Department of Electrical Engineering, National Taiwan University (NTU), Taipei, Taiwan, 2009,2016 .-{R|http://djj.ee.ntu.edu.tw/TFW.htm}- Template:Wayback
• Stéphane Mallat, A Wavelet Tour of Signal Processing Template:Wayback