Python 数据降噪处理的四种方法——均值滤波、小波变换、奇异值分解、改变binSize
github主页：https://github.com/Taot-chen
一、均值滤波
1）算法思想
?给定均值滤波窗口长度，对窗口内数据求均值，作为窗口中心点的数据的值，之后窗口向后滑动1，相邻窗口之间有重叠；边界值不做处理，即两端wid_length//2长度的数据使用原始数据。
2）Python实现
''' 均值滤波降噪：函数ava_filter用于单次计算给定窗口长度的均值滤波函数denoise用于指定次数调用ava_filter函数，进行降噪处理 '''def ava_filter(x, filt_length):N = len(x)res = []for i in range(N):if i <= filt_length // 2 or i >= N - (filt_length // 2):temp = x[i]else:sum = 0for j in range(filt_length):sum += x[i - filt_length // 2 + j]temp = sum * 1.0 / filt_lengthres.append(temp)return resdef denoise(t, x, n, filt_length):for i in range(n):res = ava_filter(x, filt_length)x = resreturn (t, res)
二、奇异值分解
1）算法思想
?任意m ? n 的矩阵A可以分解为如下形式：
A=U·sigema·V(T)
其中U、V分别是左右奇异矩阵，sigema是对角矩阵，对角线上的元素是A的奇异值从大到小的排列。
?奇异值表示的是原矩阵在其对应特征向量分量上的权重，奇异值越大，对应的特征向量在原矩阵中的权重越大。
?如果前k（k<r，r是原矩阵的秩）个奇异值数值较大，说明前k个奇异值对应的信息是原矩阵的主成分。那么可以使前k个奇异值不变，其余奇异值设置成0，再重构原矩阵，实现降噪。
2）Python实现
import numpy as np
# import random
import matplotlib.pyplot as plt
import sys
import osdef denoise(t, x):# 1、数据预处理res = int(np.sqrt(len(x)))xr = x[:res * res]delay = t[:res * res]# 2、一维数组转换为二维矩阵x2list = []for i in range(res):x2list.append(xr[i * res:i * res + res])x2array = np.array(x2list)# 3、奇异值分解U, S, V = np.linalg.svd(x2array)S_list = list(S)## 奇异值求和S_sum = sum(S)##奇异值序列归一化S_normalization_list = [x / S_sum for x in S_list]# 4、画图X = []for i in range(len(S_normalization_list)):X.append(i + 1)fig1 = plt.figure().add_subplot(111)fig1.plot(X, S_normalization_list)fig1.set_xticks(X)fig1.set_xlabel('Rank', size=15)fig1.set_ylabel('Normalize singular values', size=15)plt.show()# 5、数据重构K = 2 ## 保留的奇异值阶数for i in range(len(S_list) - K):S_list[i + K] = 0.0S_new = np.mat(np.diag(S_list))reduceNoiseMat = np.array(U * S_new * V)reduceNoiseList = []for i in range(len(x2array)):for j in range(len(x2array)):reduceNoiseList.append(reduceNoiseMat[i][j])# 6、返回结果return (delay, reduceNoiseList)
三、小波变换
1）算法思想
?将信号通过小波变换后，信号产生的小波系数含有信号的重要信息，将信号经小波分解后小波系数较大，噪声的小波系数较小，并且噪声的小波系数要小于信号的小波系数，通过选取一个合适的阀值，大于阀值的小波系数被认为是有信号产生的，应予以保留，小于阀值的则认为是噪声产生的，置为零从而达到去噪的目的。
2）Python实现
#模块调用
import numpy as np
import math
import pywt#封装成函数
def sgn(num):if (num > 0):return 1.0elif (num == 0):return 0.0else:return -1.0def wavelet_noising(new_df):data = new_dfdata = data.values.T.tolist() # 将np.ndarray()转为列表w = pywt.Wavelet('sym8')# [ca3, cd3, cd2, cd1] = pywt.wavedec(data, w, level=3) # 分解波[ca5, cd5, cd4, cd3, cd2, cd1] = pywt.wavedec(data, w, level=5) # 分解波length1 = len(cd1)length0 = len(data)Cd1 = np.array(cd1)abs_cd1 = np.abs(Cd1)median_cd1 = np.median(abs_cd1)sigma = (1.0 / 0.6745) * median_cd1lamda = sigma * math.sqrt(2.0 * math.log(float(length0), math.e))usecoeffs = []usecoeffs.append(ca5) # 向列表末尾添加对象#软硬阈值折中的方法a = 0.5for k in range(length1):if (abs(cd1[k]) >= lamda):cd1[k] = sgn(cd1[k]) * (abs(cd1[k]) - a * lamda)else:cd1[k] = 0.0length2 = len(cd2)for k in range(length2):if (abs(cd2[k]) >= lamda):cd2[k] = sgn(cd2[k]) * (abs(cd2[k]) - a * lamda)else:cd2[k] = 0.0length3 = len(cd3)for k in range(length3):if (abs(cd3[k]) >= lamda):cd3[k] = sgn(cd3[k]) * (abs(cd3[k]) - a * lamda)else:cd3[k] = 0.0length4 = len(cd4)for k in range(length4):if (abs(cd4[k]) >= lamda):cd4[k] = sgn(cd4[k]) * (abs(cd4[k]) - a * lamda)else:cd4[k] = 0.0length5 = len(cd5)for k in range(length5):if (abs(cd5[k]) >= lamda):cd5[k] = sgn(cd5[k]) * (abs(cd5[k]) - a * lamda)else:cd5[k] = 0.0usecoeffs.append(cd5)usecoeffs.append(cd4)usecoeffs.append(cd3)usecoeffs.append(cd2)usecoeffs.append(cd1)recoeffs = pywt.waverec(usecoeffs, w)return recoeffsdef denoise(x, data):data_denoising = wavelet_noising(data) #调用小波去噪函数return (x, data_denoising)
四、改变 bin size
1）算法思想
?通过改变数据的 bin size,来达到降低噪声的目的。
?改变 bin size 的时候，会导致数据长度减小，降低数据的分辨率。为了最大限度的较少原数据的有效信息的损失，在改变 bin size 的过程中，被抛弃的数据的信息也会保留在保留下来的数据中，具体实现思路是：在给定 bin size = n 的情况下，将 n 长度的数据取平均值作为该区域中心点的数据的值。之后窗口向后滑动 n ，相邻两个窗口之间不重叠。
2）Python实现
# 修改现有数据的bin：
# 即bin=3时：每三个数据，只取中间的一个数据，且这个数据的值为三个数据的平均值
# 在对纵轴进行如上处理的时候，横轴也进行相应的抽值处理：
# 第一个数据不要，第二个数据开始，每隔两个数据quyige
# 进行上述处理之前对数据进行截断处理，使数据长度为3的倍数+1，横轴数据和纵轴数据都进行截断处理
# bin=n的时候，前n//2个数据不要，后面每隔n-1个数据取一个数据，数据长度截断为n的倍数+n//2def ch_bin(x, y, bin):N = len(x)relen = N // bin * binre_x = x[:relen]re_y = y[:relen]res_x = []res_y = []i = 0while (True):if i <= bin // 2:i += 1continueelse:res_x.append(re_x[i])i += binif i >= relen - 1:breaknum = relen // binfor i in range(num):sum = 0for j in range(bin):sum += re_y[j + i * bin]res_y.append(sum * 1.0 / bin)if bin == 3:return (res_x, res_y[1:])else:return (res_x, res_y)