第二章 濾波方法介紹
2.4 離散小波轉換之介紹
小波變換的概念是由法國工程師 Morlet 在 1974 年首先提出的,它是一個時間和 頻率的局域變換,能夠有效地從訊號中提取資訊,透過伸縮或平移等運算功能,構成 𝐿2(𝑅)函數空間的一組完整的基底函數,對函數或信號進行多尺度分析。小波分析為一 種強大的數學工具,目前被廣泛地應用於信號處理的領域,可用一簡單的多分辨分析 小波樹狀圖來表示(圖 2.4)。
1988 年 Mallet 提出了多分辨分析(Multi-Resolution Analysis,MRA)的概念,將所 有的正交小波基的構造法統一起來,說明瞭小波的多解析度特性。小波分析不同於短 時傅立葉變換在單解析度上的缺點,具有多分辨分析的優點,在時域和頻域都可以處 理訊號的局部訊息。小波變換具有下列特性:時頻局部化特性、多解析度特性以及小 波基選擇的靈活性。
2.4.1 離散小波變換的推導
離散小波變換( Discrete Wavelet Transform)在數值分析和時頻分析中很有用。離散 小波轉換顧名思義就是離散的輸入以及離散的輸出。
在連續小波中,考慮小波函數:
𝛹𝑎,𝑏(𝑡) = |𝑎| − 12 𝛹(𝑡−𝑏𝑎 ) (2.21)
其中𝑏 ∈ 𝑅,𝑎 ∈ 𝑅+,且𝑎 ≠ 0。
在實際運用上,由於數據處理上的問題,連續小波必須經過離散化。限制 a >0,
所得之相容性條件變為:
10
𝐶𝛹 = ∫0∞|𝛹̂ (𝜔(𝜔̅ )̅ )|2𝑑𝜔 ̅̅̅ < ∞ (2.22)
Meyer 於 1986 年發表了具有一定衰減特性之光滑小波函數,此函數之二進制尺度 與平移足以表示𝐿2(𝑅)空間中之所有函數。將連續小波變換中的尺度因數 a 和平移因數 b 進行二進位之離散化,取 𝑎 = 𝑎0𝑗,𝑏 = 𝑘𝑎0𝑗𝑏0,其中 𝑘 , 𝑗 ∈ 𝑍。
可得離散小波函數為:
𝛹𝑗,𝑘(𝑡) = 𝑎0 − 𝑗2 𝛹(𝑡−𝑘𝑎0𝑗𝑏0
𝑎0𝑗 ) = 𝑎0 − 𝑗2 𝛹(𝑎0−𝑗− 𝑘𝑏0) (2.23) 通常我們對其進行二進制離散化,取𝑎0 =2 且 𝑏0 =1,則𝛹𝑗,𝑘(𝑡)可表示為:
𝛹𝑗,𝑘(𝑡) = 2 − 2𝑗 𝛹(2−𝑚𝑡 − 𝑛) (2.24)
此亦稱為二進制小波(Dyadic Wavelet)。
我們將離散化小波變換表示為:
𝐷𝑊𝑓(𝑗, 𝑘) = ∫ 𝑓(𝑡)𝛹𝑅 ̅𝑗,𝑘(𝑡)𝑑𝑡 (2.25)
若用內積的方式表示則寫為:
𝑫𝑾𝑓(𝑗, 𝑘) = 〈 𝒇(𝑡), 𝜳𝑗,𝑘(𝑡) 〉 𝑗, 𝑘 ∈ 𝑍 (2.26)
值得注意的是,離散小波轉換𝐷𝑊𝑗,𝑘基本上仍然將信號𝑓(𝑡)經過一系列帶通濾波器而輸 出,不同的是其帶通濾波器之中心頻率和帶寬由於 a 之離散取樣而成為一系列之離散 值,並且濾波之後的輸出也因為 b 的離散取樣而成為若干離散取樣值。
若取具正交性之小波:
〈𝛹𝑚,𝑛,𝛹𝑗,𝑘〉 = {1, 當𝑚 = 𝑗 且 𝑛 = 𝑘
0, 當𝑚 ≠ 𝑗 或 𝑛 ≠ 𝑘 , (2.27)
11
則對任一函數𝑓(𝑡) ∈ 𝐿2(𝑅),可經由小波轉換表示成:
𝑓(𝑡) = ∑𝑁𝑛=1𝐶𝑉(𝑀, 𝑛)Φ𝑀,𝑛+ ∑𝑀𝑚=1∑𝑁𝑛=1𝐶𝑊(𝑚, 𝑛)𝛹𝑚,𝑛 (2.28) 其中Φ代表尺度函數(Scaling function)。
當取不同之 m 時,即是將𝐿2(𝑅)分解成不同的空間,各個不同分解空間之關係可 以為(Mallat,1989a and 1989b):
𝐿2(𝑅) = 𝑊1⨁ 𝑉1 = 𝑊1⨁ 𝑊2⨁ 𝑉2
= 𝑊1⨁ 𝑊2⨁ ⋯ ⨁𝑊𝑀⨁ 𝑉𝑀 (2.29)
其中 𝑊𝑖⨁ 𝑉𝑖 且 𝑉𝑖 = 𝑉𝑖+1⨁𝑊𝑖+1,而 𝑉𝑖空間之訊號較𝑊𝑖低頻。𝛹𝑚,𝑛與Φ𝑚,𝑛分別代表𝑊𝑚 與 𝑉𝑚之正交基底。理論上,所選取之 M 足夠大時, 𝑉𝑀將趨近於零空間。
2.4.2 離散小波在濾波上的應用
信號和噪訊在小波域中呈現不同的狀態,兩者的小波係數會隨著尺度的不同而產 生不同的變化,隨著尺度分解的階數越多,噪聲係數的幅值會很快地衰減,而原始訊 號係數的幅值基本上維持不變。濾波的重點就在於盡可能減小由噪訊所產生的係數,
且盡可能地保留原始訊號的係數,最後將經過處理的小波係數重構回去,得到原始信 號估計的最優化。
小波濾波的三大基本假設:
(1) 噪訊經小波變換之後,大部分的小波係數趨近零。
(2) 噪訊均勻地分佈在所有係數之中。
(3) 噪訊比不至於過高[30]。
小波的信號濾波( filtering,亦可稱為去噪,de-noising)的過程可概略分為三個步 驟:
(1) 選擇一個小波分解的尺度 (階數)及小波基,進行分解計算。
(2) 選擇一個閾值(threshold),對各個分解尺度下的小波係數進行非線性處理。
(3) 將處理過後的低頻及高頻係數進行重構(小波逆變換),見圖 2.5。
普遍來說,地震訊號中的噪訊只存在於前三或四尺度的分解結果中[43],但並無特定 的判斷準則來確定小波分解的尺度,只能大概估算。
12
小波基的選擇又是另一門學問,一個良好的濾波效果,需要在緊支撐性、平滑性 以及對稱性等達到一定的水準。然而目前階段而言,這三點不太可能同時達到滿足:
緊支撐性與平滑性不可兼得,要有較好的平滑特性,則緊支撐長度需要增加;反之,
緊支撐的長度較小,雖然可以保留訊號的局部特性,但是相對的平滑性也會下降。本 文中,我們參照柳等人[43]的建議來選取最優小波基,採用 sym4 小波基。
在小波濾波的過程中,閾值對於濾波的效果具有決定性的影響,目前有許多文獻 提出了各種閾值的選取方法及理論,例如:通用閾值法是 Dohono 和 Johnstone 將多維 獨立正態變量決策理論應用於高斯噪聲模型所得出的通用閾值t = σ√2lnN[32];潘等 人驗證了基於正態分佈的自適應閾值 𝑡 = 𝑐𝜎,𝑐 ∈ [3.0,4.0][44,45]。極小化風險閾值法 則利用SURE(Stein’s Unbiased Risk Estimator)法[33,46]、交叉驗證(Cross-Validation, CV) 演算法[47]、廣義交叉驗證(Generalized Cross- Validation, GCV)演算法[48,49]對均方差 函數進行估計以得到濾波後的結果與其真實訊號之間的差異來確定閾值……等等。不 過,每個訊號的去噪效果對閾值的反應不一,目前並沒有很明確的方法可以判別何種 閾值會產生比較好的去噪效果。一般情況下,閾值的選定會採取試誤法(trial and error method)。
閾值確定後,再對小波係數進行閾值量化處理,目前主要的方法有硬閾值法及軟 閾值法,其基本概念都是去除小幅值的係數,對較大幅值的係數進行收縮或保留 [30~35]。硬閾值法因為有不連續性,通常會使得結果產生較大的方差,優點為可以保 留信號邊緣的局部特徵,缺點為易造成偽吉布斯(Gibbs)現象。軟閾值法因為收縮性質 較大,濾波結果相對平滑許多,卻容易造成濾波結果有較大的偏差,且使得邊緣模糊 等失真現象。
硬閾值法的基本處理方式為當某時刻的小波變換值大於閾值 t 時將保留原值,否 則置零,其公式為:
𝜃̂𝑗,𝑖 = {0 , |ω𝑗,𝑖| ≤ 𝑡
ω𝑗,𝑖 ,|ω𝑗,𝑖| ≤ 𝑡 , 其中 𝑡 = σ√2𝑙𝑛𝑁 (2.30)
軟閾值法則是當某時刻的小波變換值大於閾值 t 時,將其值減去 t 值,反之則將 其設為零,公式如下:
𝜃̂𝑗,𝑖 = sign(ω𝑗,𝑖) ∙ ( |ω𝑗,𝑖| − 𝑡) , 𝑡 = σ√2𝑙𝑛𝑁 (2.31)
其中sign(ω𝑗,𝑖)為符號函數,定義同式(2.20a)。
13
1995 年,Gao 提出了結合硬閾值法與軟閾值法的一種折衷形式,稱為半軟閾值法 (semisoft shrinkage) [50],1997 年並推導出基於該方法的 Minimax 閾值[51],至於上閾 值和下閾值的選取方式,請參閱文獻[50],於此不再贅述。在此法中,需要確定兩個 閾值,會增加計算上的複雜度,但可以截長補短,保留硬閾值法及軟閾值法各自的優 點。半軟閾值法的表示式如下:
𝜃̂𝑗,𝑖 = {
0 ,|ω𝑗,𝑖| ≤ 𝑡2 sign(ω𝑗,𝑖)𝑡2( |ω𝑡 𝑗,𝑖|−𝑡1)
2−𝑡1 , 𝑡1 < |ω𝑗,𝑖| ≤ 𝑡2 ω𝑗,𝑖 ,|ω𝑗,𝑖| ≤ 𝑡2
(2.32)
其中𝑡1為上閾值,𝑡2為下閾值。
14