一維信號小波去噪原理解析及python實現(xiàn)方式

更新時間：2023年06月25日 09:58:25 作者：衷科知眠

這篇文章主要介紹了一維信號小波去噪原理解析及python實現(xiàn)方式，具有很好的參考價值，希望對大家有所幫助。如有錯誤或未考慮完全的地方，望不吝賜教

一維信號小波去噪原理及python實現(xiàn)

信號去噪是經常用到的信號預處理過程，以達到在保留原有信號真是信息的基礎上盡可能低降低或者消除噪聲，獲得更高質量的信號，從而為下一步的處理奠定基礎。

去噪方法可分為時域方法與頻域方法。時域方法是指直接在原始信號上進行處理，比如均值濾波器、中值濾波器、EMD分解等方法。頻域方法是指在信號的變換域進行去噪然后再恢復到時域得到去噪后的信號，比如小波變換、傅里葉變換等方法。

無論是一維信號還是二維信號其原理都是一樣的，只不過實現(xiàn)上所有不同。這里主要記錄一下我目前需要處理的一維信號的去噪。

噪聲可分為加性噪聲和乘性噪聲兩種形式，加性噪聲一般假設其與原始信號無關，假設含噪信號為 y(i) ，其可表示為：

y(i) = x(i)+n(i)

其中 x(i) 表示原始的干凈信號， n(i) 表示噪聲信號。

而乘性噪聲則假設噪聲與信號有關，假如噪聲和信號成正比，則含噪信號 y(i) 可表示為：

y(i)=x(i)[1+n(i)]

其中第二個噪聲項是受原始信號影響的，原始信號的值越大，則噪聲越大。

加性噪聲的處理與計算都比較簡單，所以一般均假設噪聲為加性噪聲。

基于小波變換的去噪方法其具體原理這里就不闡述了，其主要思想就是噪聲在經過小波變換后的小波系數(shù)會遠遠小于原始信號的小波系數(shù)，因此可通過簡單的閾值處理就能夠去除絕大部分的噪聲，同時又不會對原始信號造成很大的損失。對具體原理感興趣的同學可以參考其他相關資料。

python 實現(xiàn)

使用python去噪，不需要我們自己動手實現(xiàn)相應的原理，安裝相應的包即可完成。

首先安裝小波變換的python包PyWavelets

pip install PyWavelets

安裝好之后就可以直接使用，我們可通過定義去噪函數(shù)將去噪過程打包，假設去噪函數(shù)名為wavelet_denoising, 其去噪過程及參數(shù)設置代碼如下：

def wavelet_denoising(signal):
    # signal: list
    wavelet_basis = 'db8'  # 此處選擇db8小波基
    w = pywt.Wavelet(wavelet_basis) 
    maxlevel = pywt.dwt_max_level(len(signal), w.dec_len) # 根據(jù)數(shù)據(jù)長度計算分解層數(shù)
    threshold = 0.1  # threshold for filtering
    coeffs = pywt.wavedec(signal, wavelet_basis, level=maxlevel)  # 將信號進行小波分解
    # 使用閾值濾除噪聲
    for i in range(1, len(coeffs)):
        coeffs[i] = pywt.threshold(coeffs[i], threshold*max(coeffs[i]))  # 將噪聲濾波
        # pywt.threshold 方法有三種閾值模式：硬閾值、軟閾值、軟硬結合
        # 默認為soft 軟閾值模式，其調用格式為
        # pywt.threshold(signal, threshold_value, mode='soft', substitute=0)
        # 處理后值為 signal/np.abs(signal) * np.maximum(np.abs(signal) - threshold_value, 0)
        # 具體可參考pywavelets官方文檔的說明
    # 重建信號
    signalrec= pywt.waverec(coeffs, wavelet_basis)
    return signalrec

通過調用wavelet_denoising方法，即可實現(xiàn)對原始信號的去噪。

值得注意的是去噪前后信號的長度會有所變化，若原始信號樣本長度為奇數(shù)，則去噪后長度加1，若原始長度為偶數(shù)，則去噪后長度不變。

通過實驗發(fā)現(xiàn)，如果要計算去噪前后信號的差，對于原始數(shù)據(jù)長度為奇數(shù)時，應該計算去噪后信號的第一個元素到倒數(shù)第二個元素這個子序列與原始信號之間的差值，也即舍去因去噪而多出的最后一個值。

如果使用去噪后的第二個元素到最后一個元素與原始信號的差，則可能出現(xiàn)偏差較大的情況。

如下圖所示，其中藍色線為去噪信號的第一個值到倒數(shù)第二個值與原始信號的差，紅色線為去噪信號的第二個值到最后一個值與原始信號的差值。

小波去噪代碼

先放代碼，代碼如下：

a = np.loadtxt('C:\\Users\\Administrator\\Desktop\\數(shù)據(jù)文件\\20201008009.txt', dtype=float)# 20201008009為一個一維列表
index = []
data = []
coffs = []
for i in range(len(a) - 1):
    X = float(i)
    Y = float(a[i])
    index.append(X)
    data.append(Y)
# create wavelet object and define parameters
w = pywt.Wavelet('db8')  # 選用Daubechies8小波
maxlev = pywt.dwt_max_level(len(data), w.dec_len)
print("maximum level is" + str(maxlev))
threshold = 0.1  # Threshold for filtering
# Decompose into wavelet components,to the level selected:
coffs = pywt.wavedec(data, 'db8', level=maxlev)  # 將信號進行小波分解
for i in range(1, len(coffs)):
    coffs[i] = pywt.threshold(coffs[i], threshold * max(coffs[i]))
datarec = pywt.waverec(coffs, 'db8')  # 將信號進行小波重構
mintime = 0
maxtime = mintime + len(data)
print(mintime, maxtime)
plt.rcParams['font.sans-serif']=['SimHei'] #使得標題輸出為中文
plt.rcParams['axes.unicode_minus'] = False #使得標題輸出為中文
plt.figure(1)
#plt.subplot(3, 1, 1)
plt.plot(index[mintime:maxtime], data[mintime:maxtime])
plt.title("原始信號")
plt.figure(2)
#plt.subplot(3, 1, 2)
plt.plot(index[mintime:maxtime], datarec[mintime:maxtime])
plt.title("小波降噪信號")
plt.figure(3)
#plt.subplot(3, 1, 3)
plt.plot(index[mintime:maxtime], data[mintime:maxtime] - datarec[mintime:maxtime])
plt.xlabel('time (s)')
plt.ylabel('error (uV)')
plt.tight_layout()
plt.show()