在之前的系列文章里，我們介紹了EEMD、CEEMD、CEEMDAN、VMD、ICEEMDAN、LMD、EWT，我們繼續(xù)補(bǔ)完該系列。

今天要講到的是小波分解，通常也就是指離散小波變換（Discrete Wavelet Transform, DWT）。在網(wǎng)上有一些介紹該方法的文章，但是總感覺不夠通俗或不夠透徹，希望讀完這篇能讓你有所收獲。

一、從小波分析到小波分解

小波分析是一種時(shí)頻域分析方法，該方法兼顧了信號(hào)在時(shí)域和頻域的信息。知乎上有一篇文章對(duì)小波分析的理解進(jìn)行了生動(dòng)的講解，建議對(duì)小波分析概念不熟的同學(xué)先看一下。咚懂咚懂咚：能不能通俗的講解下傅立葉分析和小波分析之間的關(guān)系？這篇文章中最后給出的小波變換的結(jié)果是這樣的：

圖1 連續(xù)小波變換

看起來十分厲害，不過同時(shí)會(huì)發(fā)現(xiàn)兩個(gè)問題：運(yùn)算量很大；只有數(shù)值解，沒有解析解。上述這種小波分析方法叫連續(xù)小波變換（continuous wavelet transform, CWT）。

為了減少變換運(yùn)算量，去除不必要的重復(fù)的系數(shù)，實(shí)際中使用的通常是離散小波變換（discrete wavelet transform, DWT）。

這里的“離散”指的是什么呢？

讓我們先回到小波基波（也叫母小波）的表達(dá)式：

其中s是尺度參數(shù)，表征頻率；t是位移參數(shù)，表征時(shí)間。這部分在答友的連接里也提到了。再看上一張圖，xy坐標(biāo)分別是SCALE和TRANSLATION，也就是s和t，他們在連續(xù)小波變換中是連續(xù)的。

所以，在離散小波變換中，“離散”的就是參數(shù)s和t。此時(shí)小波表達(dá)式寫為：

j和k都是整數(shù)，通常取s0=2，τ0=1。

可以看出，隨著j取值的遞增，我們可以得到一串不同的小波（子小波，也叫女兒小波...）。這些子小波的尺度參數(shù)以2的j次方的形式增長。當(dāng)使用這一系列的子小波，對(duì)一個(gè)連續(xù)函數(shù)進(jìn)行離散分析時(shí)，我們所獲得的是一組小波分析的系數(shù)，這個(gè)分析過程稱為小波系列分解。

上邊說道，尺度參數(shù)表征的是頻率，在子小波中尺度參數(shù)以2的倍數(shù)增長（即小波的“長度”被“拉長”了2倍），那么子小波對(duì)應(yīng)能檢測到的頻率值也會(huì)以1/2的倍數(shù)縮小。母小波所對(duì)應(yīng)的頻譜位于頻率譜的高端，具有最大的頻率譜范圍- 而其他的子小波的頻率譜則依次向頻譜圖的低頻端移動(dòng)，同時(shí)它們所覆蓋的頻率譜范圍也相應(yīng)地遞減。在理想的情況下，所有的濾波器應(yīng)該首尾相接互相覆蓋。

圖2 不同尺度的子小波在小波頻率譜上的覆蓋

是的，每個(gè)子小波就相當(dāng)于一個(gè)濾波器，離散小波變換的過程就是逐級(jí)濾波的過程。

具體流程是怎樣的呢？

用一句話描述就是：一組離散信號(hào)通過一系列的低通和高通濾波器，分別可以得到近似信號(hào)（用字母A表示）和細(xì)節(jié)信號(hào)（用字母D表示）。

用一張圖描述就是：

圖3 LP為低通濾波器，HP為高通濾波器，B為帶寬，2B為2倍帶寬

用一個(gè)例子來描述就是：

這段合成信號(hào)是由多個(gè)成分疊加而成的：一個(gè)10 Hz的低頻正弦波和一個(gè)振幅為0.5的50 Hz高頻正弦波，這兩者代表信號(hào)的周期性成分。同時(shí)，信號(hào)中還包括了一個(gè)線性趨勢項(xiàng)，表示信號(hào)的非周期性變化。除此之外，信號(hào)還疊加了標(biāo)準(zhǔn)差為0.5的高斯白噪聲，為信號(hào)添加了隨機(jī)性。這里使用了 dB4（第4級(jí)Daubechies小波）作為去除噪音操作的母小波。

原始信號(hào)

一階小波分解的結(jié)果為：

一階小波分解的近似信號(hào)（低通結(jié)果）

一階小波分解的細(xì)節(jié)信號(hào)（高通結(jié)果）

二階小波分解的結(jié)果為（即對(duì)A1信號(hào)做分解）：

二階小波分解的近似信號(hào)（低通結(jié)果）

二階小波分解的細(xì)節(jié)信號(hào)（高通結(jié)果）

三階小波分解的結(jié)果為（即對(duì)A2信號(hào)做分解）：

三階小波分解的近似信號(hào)（低通結(jié)果）

三階小波分解的細(xì)節(jié)信號(hào)（高通結(jié)果）

四階小波分解的結(jié)果為（即對(duì)A3信號(hào)做分解）：

四階小波分解的近似信號(hào)（低通結(jié)果）

四階小波分解的近似信號(hào)（低通結(jié)果）

至此我們已經(jīng)能夠得到較好的低通濾波結(jié)果了?？梢钥吹皆夹盘?hào)被逐級(jí)的，無遺漏地進(jìn)行了高、低通濾波，且越接近低頻分段越細(xì)，幾乎想要哪個(gè)頻段的特征都能得到，因而這個(gè)方法有個(gè)霸氣的名字，叫filter banks。

小波分解的多尺度可以類比為我們使用不同的“放大鏡”去觀察一個(gè)物體。想象一下你手里有一張非常復(fù)雜的畫，畫面上有大的物體，如山脈、樹木，但也有非常細(xì)小的細(xì)節(jié)，如葉子上的紋理或昆蟲的觸角。
粗尺度（低分辨率）：當(dāng)你使用低倍的放大鏡（或者站得很遠(yuǎn)）去看這幅畫時(shí)，你可以看到大的物體，如山脈和樹木，但可能看不到細(xì)小的紋理或昆蟲。 在小波分解中，這就像我們查看信號(hào)的低頻部分，捕獲其主要的、寬泛的特征。
細(xì)尺度（高分辨率）：現(xiàn)在，如果你換一個(gè)高倍的放大鏡（或者走近一些）去看同一幅畫，你可能會(huì)失去對(duì)整體的感知，但可以清晰地看到葉子上的紋理或昆蟲的觸角等細(xì)節(jié)。 在小波分解中，這就像我們查看信號(hào)的高頻部分，捕獲其細(xì)節(jié)和快速的變化。
小波分解的美妙之處在于，它同時(shí)提供了多個(gè)尺度的視角，讓我們既可以看到信號(hào)的整體特征，又可以看到其細(xì)節(jié)。這就像我們可以同時(shí)擁有多個(gè)不同倍率的放大鏡，讓我們在需要的時(shí)候選擇合適的一個(gè)來觀察畫面。

二、小波分解更深一步的理解

我們先看一下下邊這張圖：

小波分解中的c與l

乍一看這張圖與圖3比較相像，不過仔細(xì)看下邊的兩個(gè)方框，即左側(cè)分別標(biāo)識(shí)了c和l的位置。

（一）關(guān)于小波分解系數(shù)

c的那一行是指的就是小波分解的向量。需要注意這個(gè)向量并不是上邊幾張圖里的近似信號(hào)和細(xì)節(jié)信號(hào)，這里邊存儲(chǔ)的是小波分解的系數(shù)。小波系數(shù)是沒有量綱單位的結(jié)果，需要經(jīng)過重構(gòu)這些系數(shù)得到實(shí)際有量綱的信號(hào)。

在我們之前講過的小波閾值去噪方法就是針對(duì)這些小波系數(shù)，c這行的每個(gè)獨(dú)立的方框（比如cA3）都可以重構(gòu)到時(shí)域成為對(duì)應(yīng)的濾波后的信號(hào)，也可以幾個(gè)方框共同重構(gòu)，這就是小波分解分量的篩選重構(gòu)過程，這個(gè)過程是有實(shí)際工程應(yīng)用意義的，比如可以實(shí)現(xiàn)濾波。

（二）小波分解和“其他類EMD分解方法”的區(qū)別

這里指的“其他類EMD方法”包括了EMD、EEMD、CEEMD、CEEMDAN、ICEEMDAN等。

小波分解與EMD分解最大的不同是源于分解機(jī)理的。

小波分解的分解結(jié)構(gòu)是有包含關(guān)系的，每一層級(jí)的近似信號(hào)都要再分解為下一級(jí)的近似信號(hào)和細(xì)節(jié)信號(hào)，（這也直接導(dǎo)致在信號(hào)重構(gòu)時(shí)不是那么隨意），下圖是一個(gè)典型信號(hào)的分解結(jié)果，我用紅色箭頭標(biāo)出了其包含關(guān)系。

主要低頻信號(hào)出現(xiàn)在每一個(gè)近似信號(hào)里

EMD通過連續(xù)地提取信號(hào)的局部極值，然后求解其上下包絡(luò)，進(jìn)而得到IMFs。這些IMFs是并列的關(guān)系，可以直接對(duì)選定的分量進(jìn)行相加來實(shí)現(xiàn)重構(gòu)，這種操作就比較直觀了。

每個(gè)IMF分量都是由原始信號(hào)直接分解而來

那么這種分解機(jī)理的區(qū)別，會(huì)帶來分解效果和應(yīng)用上的什么不同呢？

我認(rèn)為有以下幾點(diǎn)：

1. 上邊講到小波分解方法是多尺度的，高頻分量在多個(gè)層級(jí)上被逐步細(xì)化剝離開來，讓我們可以用更精準(zhǔn)的手術(shù)刀切割出特征信號(hào)段，這是EMD所不具備的特點(diǎn)。
2. EMD可能會(huì)受到模態(tài)混疊的影響，導(dǎo)致不同模式的信號(hào)成分被混在一起；小波分解方法不同層級(jí)覆蓋的頻率范圍不同，模態(tài)混疊現(xiàn)象會(huì)大大減少。
3. EMD的分解階數(shù)是自適應(yīng)的，換句話說無法人為干預(yù)；小波分解的分解層數(shù)是可以指定的，這方面靈活性更強(qiáng)。

三、小波分解DWT的MATLAB代碼實(shí)現(xiàn)

小波分解的代碼在網(wǎng)上可以找到一些，但是用起來不太趁手。

按照“類EMD”系列的代碼的統(tǒng)一風(fēng)格，筆者進(jìn)行了封裝，封裝后的函數(shù)有三個(gè)，分別用于實(shí)現(xiàn)繪制小波分解圖、小波分解各分量及頻譜對(duì)應(yīng)圖，以及重構(gòu)信號(hào)并繪制信號(hào)重構(gòu)圖。

（一）生成仿真信號(hào)

%% 1.生成仿真信號(hào)
Fs = 1000; % 采樣頻率
t = 0:1/Fs:1-1/Fs; % 時(shí)間向量
% 創(chuàng)建一個(gè)合成信號(hào)：包含不同頻率的正弦波、趨勢和噪聲
signal = cos(2*pi*10*t) + 0.5*sin(2*pi*50*t) + t + 0.5*randn(size(t));
figure('color','white')
plot(t,signal,'k');xlabel('時(shí)間');ylabel('幅值') %繪制原始信號(hào)

待分解的原始信號(hào)

（二）小波分解圖

%% 2.繪制DWT分解圖
waveletType = 'db4';    %小波名稱,可選范圍參考這里：https://ww2.mathworks.cn/help/wavelet/ref/wfilters.html?searchHighlight=wname&s_tid=srchtitle_wname_2#d123e130597
decompositionLevel = 4; %小波分解水平，正整數(shù)
[a,d] = pDWT(signal, decompositionLevel, waveletType); % 調(diào)用函數(shù)進(jìn)行分解和畫圖

只需要設(shè)置小波名稱和小波分解水平，然后調(diào)用pDWT函數(shù)即可（函數(shù)獲取方法見文末）

此時(shí)可以畫出如下圖：

（三）小波分解及頻譜圖

%% 3.繪制DWT分解圖及頻譜圖
waveletType = 'db4';    %小波名稱,可選范圍參考這里：https://ww2.mathworks.cn/help/wavelet/ref/wfilters.html?searchHighlight=wname&s_tid=srchtitle_wname_2#d123e130597
decompositionLevel = 4; %小波分解水平，正整數(shù)
[a,d] = pDWTandFFT(signal, decompositionLevel, waveletType); % 調(diào)用函數(shù)進(jìn)行分解和畫圖（及頻譜圖）

只需要設(shè)置小波名稱和小波分解水平，然后調(diào)用pDWTandFFT函數(shù)即可（函數(shù)獲取方法見文末）

此時(shí)可以畫出如下圖：

（四）小波分解重構(gòu)及畫圖

%% 4.重構(gòu)信號(hào)并繪制DWT重構(gòu)圖
waveletType = 'db4';    %小波名稱,可選范圍參考這里：https://ww2.mathworks.cn/help/wavelet/ref/wfilters.html?searchHighlight=wname&s_tid=srchtitle_wname_2#d123e130597
decompositionLevel = 4; %小波分解水平，正整數(shù)
approxLevels = 4;%所選的近似分量
detailLevels = [3 4];%所選的細(xì)節(jié)分量
combined_signal = rDWT(signal, waveletType, decompositionLevel, approxLevels, detailLevels);

為了應(yīng)用小波分解結(jié)果（比如濾波），很多時(shí)候要對(duì)分解結(jié)果重構(gòu)。

重構(gòu)的操作相對(duì)復(fù)雜一些，一來這部分代碼一不小心可能就會(huì)寫錯(cuò)；二來重構(gòu)選擇分量的時(shí)候也容易出錯(cuò)。

對(duì)于第一個(gè)問題，我寫了一個(gè)重構(gòu)的封裝函數(shù)（就是上邊這段演示的，封裝函數(shù)為rDWT），只需要選擇想要重構(gòu)的近似分量和細(xì)節(jié)分量就行。

對(duì)于第二個(gè)問題，大家只需要注意兩個(gè)問題：1.重構(gòu)選擇近似分量的時(shí)候，近似分量approxLevels只能選擇一個(gè)層級(jí)，比如可以讓approxLevels=3或者=4，但是不能=[3,4]，否則低頻分量就被重構(gòu)了兩次，重構(gòu)后的數(shù)據(jù)會(huì)超出原始數(shù)據(jù)大小；2.選擇了高層級(jí)的近似分量后，就不能在選擇低于他的層級(jí)的細(xì)節(jié)分量，比如如果設(shè)置了approxLevels=3，就不能再選擇讓detailLevels =4。

調(diào)用上述分解重構(gòu)及畫圖函數(shù)，需要設(shè)置小波名稱和小波分解水平，然后調(diào)用pDWTandFFT函數(shù)即可（函數(shù)獲取方法見文末）

此時(shí)可以畫出如下圖：

上邊提到了三個(gè)封裝函數(shù)，分別是：

• pDWT（小波分解并畫圖的程序）
• pDWTandFFT（繪制DWT分解圖及頻譜圖）
• rDWT（重構(gòu)信號(hào)并繪制DWT重構(gòu)圖）

上邊的測試代碼和封裝函數(shù)，包括工具箱都可以在公眾號(hào)khscience（看海的城堡）中回復(fù)"DWT"獲取，EMD、EEMD、CEEMD、CEEMDAN、ICEEMDAN、VMD以及HHT相關(guān)的程序也有，編程不易，感謝支持~關(guān)于EMD、EEMD、CEEMD、VMD和HHT的相關(guān)介紹可以看這里：

Mr.看海：這篇文章能讓你明白經(jīng)驗(yàn)?zāi)B(tài)分解（EMD）——EMD在MATLAB中的實(shí)現(xiàn)方法

Mr.看海：希爾伯特譜、邊際譜、包絡(luò)譜、瞬時(shí)頻率/幅值/相位——Hilbert分析衍生方法及MATLAB實(shí)現(xiàn)

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第一篇）——EEMD

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第二篇）——CEEMD

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第三篇）——CEEMDAN

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第四篇）——VMD

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第五篇）——ICEEMDAN

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第六篇）——LMD

Mr.看海：類EMD的“信號(hào)分解方法”及MATLAB實(shí)現(xiàn)（第七篇）——EWT