欢迎您访问程序员文章站本站旨在为大家提供分享程序员计算机编程知识!
您现在的位置是: 首页  >  IT编程

荐 python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)

程序员文章站 2022-03-30 10:16:46
python 数据、曲线平滑处理——方法总结Savitzky-Golay 滤波器实现曲线平滑插值法对折线进行平滑曲线处理基于Numpy.convolve实现滑动平均滤波数据平滑处理——log()和exp()函数问题描述:在寻找曲线的波峰、波谷时,由于数据帧数多的原因,导致生成的曲线图噪声很大,不易寻找规律。如下图:由于高频某些点的波动导致高频曲线非常难看,为了降低噪声干扰,需要对曲线做平滑处理,让曲线过渡更平滑。常见的对曲线进行平滑处理的方法包括: Savitzky-Golay 滤波器、插值法等。...


问题描述:

在寻找曲线的波峰、波谷时,由于数据帧数多的原因,导致生成的曲线图噪声很大,不易寻找规律。如下图:

荐
                                                        python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)

由于高频某些点的波动导致高频曲线非常难看,为了降低噪声干扰,需要对曲线做平滑处理,让曲线过渡更平滑。常见的对曲线进行平滑处理的方法包括: Savitzky-Golay 滤波器、插值法等。

Savitzky-Golay 滤波器实现曲线平滑

对曲线进行平滑处理,通过Savitzky-Golay 滤波器,可以在scipy库里直接调用,不需要再定义函数。

代码语法:

python中Savitzky-Golay滤波器调用如下:

y_smooth = scipy.signal.savgol_filter(y,53,3)  
# 亦或
y_smooth2 = savgol_filter(y, 99, 1, mode= 'nearest')

# 备注:
y:代表曲线点坐标(x,y)中的y值数组
window_length:窗口长度,该值需为正奇整数。例如:此处取值53
k值:polyorder为对窗口内的数据点进行k阶多项式拟合,k的值需要小于window_length。例如:此处取值3
mode:确定了要应用滤波器的填充信号的扩展类型。(This determines the type of extension to use for the padded signal to which the filter is applied.

调参规律:

现在看一下window_length和k这两个值对曲线的影响。

(1)window_length对曲线的平滑作用: window_length的值越小,曲线越贴近真实曲线;window_length值越大,平滑效果越厉害(备注:该值必须为正奇整数)。

(2)k值对曲线的平滑作用: k值越大,曲线越贴近真实曲线;k值越小,曲线平滑越厉害。另外,当k值较大时,受窗口长度限制,拟合会出现问题,高频曲线会变成直线。

典型范例:

# 用于生成问题描述中示例曲线的代码如下:
import numpy as np

Size = 100
x = np.linspace(1, Size,Size)

data = np.random.randint(1, Size, Size)
print(data)
>>>
array([33, 19,  2, 56, 11, 38,  3, 45,  4, 20, 59,  3, 71, 94, 85, 23, 93,
       65, 23, 99, 38, 43, 41, 96, 46, 56, 79,  1, 38, 90, 97, 88,  2, 72,
       25, 51, 70, 42,  4, 86, 26, 44, 40, 49,  5, 37, 10, 99, 40, 88, 34,
       99, 26,  6, 20, 37, 22, 88, 13, 68, 54, 95, 15,  4, 58, 20, 51, 89,
       81,  1, 41, 21, 17, 52, 84, 46, 76, 44, 90, 72, 96, 32, 12, 50, 81,
       64, 67, 99, 42, 35, 70, 79, 21, 51, 56, 10, 23, 21, 25, 64])
      
# 可视化图线
plt.plot(x, data)

# 使用Savitzky-Golay 滤波器后得到平滑图线
from scipy.signal import savgol_filter

y = savgol_filter(data, 5, 3, mode= 'nearest')
# 可视化图线
plt.plot(x, y, 'b', label = 'savgol')

下图所示为经过Savitzky-Golay 滤波器处理的示例曲线效果图:
荐
                                                        python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)
原理剖析:

scipy官方帮助文档里可以看到关于Savitzky-Golay 滤波器中关于savgol_filter()函数的详细说明。

以下是关于Savitzky-Golay平滑滤波的简单介绍(引至博客:Python 生成曲线进行快速平滑处理):

Savitzky-Golay平滑滤波是光谱预处理中的常用滤波方法,其 核心思想:是对一定长度窗口内的数据点进行k阶多项式拟合,从而得到拟合后的结果。 对它进行离散化处理后,S-G 滤波其实是一种移动窗口的加权平均算法,但是其加权系数不是简单的常数窗口,而是通过在滑动窗口内对给定高阶多项式的最小二乘拟合得出。

Savitzky-Golay平滑滤波被广泛地运用于数据流平滑除噪,是一种在时域内基于局域多项式最小二乘法拟合的滤波方法。这种滤波器的 最大特点:在滤除噪声的同时可以确保信号的形状、宽度不变。

使用平滑滤波器对信号滤波时,实际上是拟合了信号中的低频成分,而将高频成分平滑出去了。 如果噪声在高频端,那么滤波的结果就是去除了噪声,反之,若噪声在低频段,那么滤波的结果就是留下了噪声。

总之,平滑滤波是光谱分析中常用的预处理方法之一。用Savitzky-Golay方法进行平滑滤波,可以提高光谱的平滑性,并降低噪音的干扰。S-G平滑滤波的效果,随着选取窗宽不同而不同,可以满足多种不同场合的需求。

参考链接:Savitzky-Golay平滑滤波的python实现

插值法对折线进行平滑曲线处理

实现所需的库: numpy、scipy、matplotlib

插值法的常见实现方法:

  • nearest:最邻近插值法
  • zero:阶梯插值
  • slinear:线性插值
  • quadraticcubic:2、3阶B样条曲线插值

拟合和插值的区别:

1、插值:简单来说,插值就是根据原有数据进行填充,最后生成的曲线一定过原有点。

2拟合:拟合是通过原有数据,调整曲线系数,使得曲线与已知点集的差别(最小二乘)最小,最后生成的曲线不一定经过原有点。

代码语法:

通过执行from scipy.interpolate import make_interp_spline,导入make_interp_spline模块,之后调用make_interp_spline(x, y)(x_smooth)函数实现。

官方帮助文档:scipy.interpolate.make_interp_spline

典型范例:

import numpy as np
from matplotlib import pyplot as plt
from scipy.interpolate import make_interp_spline

x = np.array([6, 7, 8, 9, 10, 11, 12])
y = np.array([1.53E+03, 5.92E+02, 2.04E+02, 7.24E+01, 2.72E+01, 1.10E+01, 4.70E+00])
x_smooth = np.linspace(x.min(), x.max(), 300)  # np.linspace 等差数列,从x.min()到x.max()生成300个数,便于后续插值
y_smooth = make_interp_spline(x, y)(x_smooth)
plt.plot(x_smooth, y_smooth)
plt.show()

上述执行代码的效果如下:
荐
                                                        python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)
参考链接:python利用插值法对折线进行平滑曲线处理

基于Numpy.convolve实现滑动平均滤波

滑动平均概念:

滑动平均滤波法 (又称:递推平均滤波法),它把连续取N个采样值看成一个队列 ,队列的长度固定为N ,每次采样到一个新数据放入队尾,并扔掉原来队首的一次数据(先进先出原则) 。把队列中的N个数据进行算术平均运算,就可获得新的滤波结果。

N值的选取:流量,N=12;压力:N=4;液面,N=4 ~ 12;温度,N=1~4

滑动平均的优缺点:

优点: 对周期性干扰有良好的抑制作用,平滑度高,适用于高频振荡的系统。

缺点: 灵敏度低,对偶然出现的脉冲性干扰的抑制作用较差,不易消除由于脉冲干扰所引起的采样值偏差,不适用于脉冲干扰比较严重的场合,比较浪费RAM 。

滑动平均的数学原理:

滑动平均滤波法计算类似一维卷积的工作原理,滑动平均的N就对应一维卷积核大小(长度)。 不过区别在于:

(1)步长会有些区别,滑动平均滤波法滑动步长为1,而一维卷积步长可以自定义;
(2)一维卷积的核参数是需要更新迭代的,而滑动平均滤波法核参数都是1。

我们应该怎么利用这个相似性呢?其实也很简单,只需要把一维卷积核大小(长度)和N相等,步长设置为1,核参数都初始为1就可以了。由于一维卷积计算速度快,因此我们可以使用一维卷积来快速高效地实现这个功能。

代码语法:

通过Numpy库中的convolve()函数可以实现这些功能。

def np_move_avg(a,n,mode="same"):
    return(np.convolve(a, np.ones((n,))/n, mode=mode))

Numpy.convolve函数:numpy.convolve函数官方文档

numpy.convolve(a, v, mode=‘full’)  # 返回a和v的离散线性卷积。

参数说明:

  • a:(N,)输入的第一个一维数组
  • v:(M,)输入的第二个一维数组
  • mode:{‘full’, ‘valid’, ‘same’}参数可选,该参数指定np.convolve函数如何处理边缘。

mode可能的三种取值情况:
full’ 默认值,返回每一个卷积值,长度是N+M-1,在卷积的边缘处,信号不重叠,存在边际效应。
‘same’ 返回的数组长度为max(M, N),边际效应依旧存在。
‘valid’  返回的数组长度为max(M,N)-min(M,N)+1,此时返回的是完全重叠的点。边缘的点无效。

和一维卷积参数类似,a就是被卷积数据,v是卷积核大小。

原理剖析:

滑动平均值是卷积数学运算的一个例子。对于滑动平均值,沿着输入滑动窗口并计算窗口内容的平均值。对于离散的1D信号,卷积是相同的,除了代替计算任意线性组合的平均值,即将每个元素乘以相应的系数并将结果相加。那些系数,一个用于窗口中的每个位置,有时称为卷积核。现在,N值的算术平均值是(x1+x2+...+xN)/N(x_1 + x_2 + ... + x_N) / N,所以相应的内核是(1/N,1/N,...,1/N)(1/N, 1/N, ..., 1/N),这正是我们通过使用得到的np.ones((N,))/N

np.convolve函数中通过mode参数指定如何处理边缘。下面是一个说明模式不同取值之间差异的图:

import numpy as np
import matplotlib.pyplot as plt
 
def np_move_avg(a,n,mode="same"):
    return(np.convolve(a, np.ones((n,))/n, mode=mode))
 
modes = ['full', 'same', 'valid']

for m in modes:
    plt.plot(np_move_avg(np.ones((200,)), 50, mode=m))
 

plt.axis([-10, 251, -.1, 1.1])
 
plt.legend(modes, loc='lower center')
 
plt.show()

荐
                                                        python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)
参考链接:
[开发技巧]·Python极简实现滑动平均滤波(基于Numpy.convolve)
numpy中的convolve的理解

典型范例:

# 实现数据可视化中的数据平滑
import numpy as np
import matplotlib.pylab as plt
 
'''
其它的一些知识点:
raise:当程序发生错误,python将自动引发异常,也可以通过raise显示的引发异常
一旦执行了raise语句,raise语句后面的语句将不能执行
'''
 
def moving_average(interval, windowsize):
    window = np.ones(int(windowsize)) / float(windowsize)
    re = np.convolve(interval, window, 'same')
    return re
 
def LabberRing():
    t = np.linspace(-4, 4, 100)   # np.linspace 等差数列,从-4到4生成100个数
    print('t=', t)
 # np.random.randn 标准正态分布的随机数,np.random.rand 随机样本数值
    y = np.sin(t) + np.random.randn(len(t)) * 0.1   # 标准正态分布中返回1个,或者多个样本值
    print('y=', y)
    
    plt.plot(t, y, 'k')     # plot(横坐标,纵坐标, 颜色)
    
    y_av = moving_average(y, 10)
    plt.plot(t, y_av, 'b')
    plt.xlabel('Time')
    plt.ylabel('Value')
    # plt.grid()网格线设置
    plt.grid(True)
    plt.show()
    return
 
LabberRing()  # 调用函数

黑色曲线为原始折线,蓝色曲线为通过np.convolve函数进行平滑处理后的曲线效果。
荐
                                                        python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)
参考链接:数据可视化中的数据平滑

数据平滑处理——log()和exp()函数

log()函数的使用就像是将一个数据压缩到了一个区间,与数据的标准化类似

在数据预处理时,引入log()函数的优点:

1、在数据预处理时首先可以对偏度比较大的数据用log()函数进行转化,使其更加服从高斯分布,此步处理可能会使我们后续的分类结果得到一个更好的结果。

2、在数据、曲线平滑处理时,同样需要使用log()函数,从而避免复值的问题 ,否则,将导致模型的结果总是达不到一定的标准。(备注: 复值指一个自变量对应多个因变量

下面再说说它的 逆运算 exp()函数

由于前面使用过log()将数据进行了压缩,所以最后需要记得将预测出的平滑数据进行一个还原,而还原过程就是log()函数的逆运算exp()函数

常用的数学表达式:
log1p=log(x+1)ln(x+1)log1p = log(x+1)即:ln(x+1)
expm1=exp(x)1expm1 = exp(x)-1

log1p函数有它存在的意义,即保证了x数据的有效性,当x很小时(如 两个数值相减后得到x=1016x = 10^{-16}),由于太小超过数值有效性,用log(x+1)log(x+1)计算得到结果为0,换作log1p则计算得到一个很小却不为0的结果,这便是它的意义(好像是用泰勒公式来展开运算的,不确定)。
同样的道理对于expm1,当x特别小,就会急剧下降出现如上问题,甚至出现错误值。
在最开始看到这样的处理方式的时候,不是很理解包括为什么是逆运算(一下子没有想到),后来慢慢摸索就优点清晰了,比如为什么两这是逆运算(简单处理):
logxlogx是e为底的对数,exe^x是e为底的指数,根据对数的规则,再进行变换推导可以得到:elogex=xe^{log_e^x}=x
可以看到x经过对数的处理后,再经过指数处理再次得到x,这里对两者的逆运算做了简单的介绍。
荐
                                                        python 数据、曲线平滑处理——方法总结(Savitzky-Golay 滤波器、make_interp_spline插值法和convolve滑动平均滤波)
另外RMSLE(均方根对数误差)会更多的惩罚欠拟合,所以在使用该误差定义时我们也可以用到上面的函数:
np.loglp计算加一后的对数,其逆运算是np.expm1;
采用此误差函数时,可以先对原始数据做np.log1p,再使用RMSE。

参考链接:数据平滑处理——log1p()和exmp1()

本文地址:https://blog.csdn.net/weixin_42782150/article/details/107176500