离散小波变换

在數值分析和泛函分析領域中，**離散小波轉換**（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]$

　
: 可繼續延伸

 

　
: : $\vdots $

　
: $\vdots $

架構中的第 $\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[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 IDFT_{N+L-1}\left[DFT_{N+L-1}(x[n])DFT_{N+L-1}(y[n])\right]$

其中，

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

 $\ DFT_{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{aligned}{\frac {3}{2}}S(N_{1}+L-1)\log _{2}{(N_{1}+L-1)}&\approx {\frac {3}{2}}SN_{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{aligned}}
$$

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

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

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

$$
{\begin{aligned}Complexity&\approx \left(N+{\frac {N}{2}}+{\frac {N}{4}}+{\frac {N}{8}}+...+2\right){\frac {3}{2}}\log _{2}N_{1}\\&=(2N-2){\frac {3}{2}}\log _{2}N_{1}\\&\approx 3N\log _{2}N_{1}\\\end{aligned}}
$$

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

$$
{\begin{aligned}Complexity&\approx \left(N+2{\frac {N}{2}}+4{\frac {N}{4}}+8{\frac {N}{8}}+...+{\frac {N}{2}}2\right){\frac {3}{2}}\log _{2}N_{1}\\&=(N\log _{2}N){\frac {3}{2}}\log _{2}N_{1}\\\end{aligned}}
$$

(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{aligned}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{aligned}}
$$

(5) 二維離散小波轉換之複雜度，使用 「分段摺積 (sectioned convolution)」的技巧：

假設原始尺寸為 $M\times N$ ，則每一小部分的尺寸為 $M_{1}\times N_{1}$

$$
{\begin{aligned}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{aligned}}
$$

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

(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{aligned}totalcomplexity&=\left(MN+{\frac {MN}{4}}+{\frac {MN}{16}}+...\right){\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{aligned}}
$$

2. 若 $\ x_{a,H_{1}}[n],x_{a,H_{2}}[n],x_{a,H_{3}}[n]$ 也細分時，總複雜度為：

$$
{\begin{aligned}totalcomplexity&=\left(MN+4{\frac {M}{2}}{\frac {N}{2}}+16{\frac {M}{4}}{\frac {N}{4}}+...\right){\frac {3}{2}}\log _{2}M_{1}N_{1}\\&=\left[MN\log _{2}(min(M,N))\right]{\frac {3}{2}}\log _{2}M_{1}N_{1}\\\end{aligned}}
$$

## 重建 (Reconstruction)

使用離散小波轉換，將訊號個別通過一個低通濾波器和一個高通濾波器，得到訊號的高低頻成分，而在重建 (Reconstruction) 原始訊號的過程，也就是離散小波的逆轉換 (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)={\tfrac {1}{2}}[X(z)G(z)+X(-z)G(-z)]G_{1}(z)+{\tfrac {1}{2}}[X(z)H(z)+X(-z)H(-z)]H_{1}(z)$
```

```
$={\tfrac {1}{2}}[G(z)G_{1}(z)+H(z)H_{1}(z)]X(z)+{\tfrac {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$
```

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

```
${\binom {G_{1}(z)}{H_{1}(z)}}={\frac {2}{det(H_{m}(z))}}{\binom {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. 滿足完整重建要條件， ${\binom {G_{1}(z)}{H_{1}(z)}}={\frac {2}{det(H_{m}(z))}}{\binom {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)$

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

 * $G_{1}(z)=G(z)z^{-k}$

$$
\Longrightarrow $ $g_{1}(n)=g(n-k)
$$

 * $H_{1}(z)=-G(-z)z^{-k}$

$$
\Longrightarrow $ $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})$

$$
\Longrightarrow $ $g(n)=g_{1}(-n)
$$

 * $H(z)=-z^{k}G_{1}(-z)$

$$
\Longrightarrow $ $h(n)=(-1)^{n}g_{1}(n+k)
$$

 * $H_{1}(z)=-z^{k}G_{1}(-z^{-1})$

$$
\Longrightarrow $ $h_{1}(n)=(-1)^{n}g_{1}(-n+k)
$$

 $det(H_{m}(z))=G(z)H(-z)-H(z)G(-z)=2z^{k}$ ，且 $k$ 為奇數。

多數的 wavelet 屬於 orthonormal wavelet。

## 其他應用

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

* 邊緣偵測：使用高通濾波器，將小波的低頻濾掉，即保留 $x_{1,HL}[m,n]$ 或 $x_{1,LH}[m,n]$ 而捨棄其他部分。
* R 語言小波分析 wavelet （頁面存檔備份，存於網際網路檔案館）
* 作為 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 .http://djj.ee.ntu.edu.tw/TFW.htm （頁面存檔備份，存於網際網路檔案館）
* Stéphane Mallat, A Wavelet Tour of Signal Processing （頁面存檔備份，存於網際網路檔案館）