二维数据平滑,滤波方法比较:移动平均平滑、中值滤波、指数加权移动平均、Savitzky-Golay、高斯滤波、卡尔曼滤波


前言

数据平滑是一种数据处理技术,主要目的是消除数据中的噪声或不规则波动,使数据更加平滑和易于分析。由于获取的原始数据或者计算后的数据,由于环境、设备等影响,数据会产生随机波动、噪声或异常值等情况,所以需要通过数据平滑对数据进行处理,使其变得更加平缓、连续或更容易被处理。


一、移动平均平滑

移动平均平滑是基于平均值的概念,通过计算序列中每个数据点周围的一定数量的数据点的平均值,来平滑时间序列中的噪声和波动,从而更清晰地观察序列的趋势和周期性。

计算方法
移动平均就是在数据中,按顺序依次选定一个固定长度的区间(窗口),然后计算这个窗口内数据的平均值,以平均值代替当前值。随着窗口不断向前滑动,每次都计算出新窗口内数据的平均值。

实现代码,我这边读取的是本地csv文件中的数据:

import pandas as pd
import matplotlib.pyplot as plt  

# 读取数据
csv_file1 = "outputdata.csv"
'''
csv文件中的数据格式为:
x,y
6.15,1.87
6.13,1.54
6.05,1.46
...

'''
df1 = pd.read_csv(csv_file1)

window_size = 5  # 窗口大小
df1['x_smoothed_ma'] = df1['x'].rolling(window=window_size, center=True).mean()   # x
df1['y_smoothed_ma'] = df1['y'].rolling(window=window_size, center=True).mean()   # y

plt.figure(figsize=(10, 6))  # 设置图形大小  
plt.plot(df1['x'], df1['y'], linewidth=4, label='original')  # 绘制原始数据  
plt.plot(df1['x_smoothed_ma'], df1['y_smoothed_ma'], linewidth=2, label='ma')  # 绘制滤波数据 

# 添加图例  
plt.legend()  
  
# 设置纵横比相等,这取决于数据的具体范围  
plt.gca().set_aspect('equal', adjustable='box')

# 添加标题和坐标轴标签  
plt.title('Data Smoothing')  
plt.xlabel('X Axis')  
plt.ylabel('Y Axis')  
plt.grid(True)

# 显示图形  
plt.show()

平滑效果:
在这里插入图片描述

二、中值滤波

中值滤波和移动平均平滑算法一致,只不过是将窗口内的均值换成中值,其他步骤一致。

实现代码:

import pandas as pd
import matplotlib.pyplot as plt  

# 读取数据
csv_file1 = "outputdata.csv"
'''
csv文件中的数据格式为:
x,y
6.15,1.87
6.13,1.54
6.05,1.46
...

'''
df1 = pd.read_csv(csv_file1)

window_size = 5  # 窗口大小
# 中值滤波
df1['x_smoothed_median'] = df1['x'].rolling(window=window_size, center=True).median()
df1['y_smoothed_median'] = df1['y'].rolling(window=window_size, center=True).median()

plt.figure(figsize=(10, 6))  # 设置图形大小  
plt.plot(df1['x'], df1['y'], linewidth=4, label='original')  # 绘制原始数据  
plt.plot(df1['x_smoothed_median'], df1['y_smoothed_median'], linewidth=2, label='ma')  # 绘制滤波数据 

# 添加图例  
plt.legend()  
  
# 设置纵横比相等,这取决于数据的具体范围  
plt.gca().set_aspect('equal', adjustable='box')

# 添加标题和坐标轴标签  
plt.title('Data Smoothing')  
plt.xlabel('X Axis')  
plt.ylabel('Y Axis')  
plt.grid(True)

# 显示图形  
plt.show()

平滑效果:
在这里插入图片描述
在移动平均平滑和中值滤波过程中在首尾会出现空值,主要有以下原因:

对于移动平均平滑,当设置了 center=True 时,在序列的开头和结尾部分,由于没有足够的数据点来进行完整窗口的计算,所以会产生空值。例如,对于窗口大小为 5 的移动平均,在序列的前 2 个点和后 2 个点就无法计算出中心的移动平均值,从而导致这些位置出现空值。

对于中值滤波,同样当设置了 center=True 且窗口大小较大时,在序列的开头和结尾部分也会因为数据不足而无法计算中值,进而产生空值。

如果不希望出现空值,可以考虑以下几种方法:

①不使用 center=True 参数,但这样可能会导致平滑后的数据在边缘处有一定的偏差。

②在进行平滑处理之前,对数据的开头和结尾进行适当的填充,例如使用固定值或其他合理的填充方法,然后再进行平滑操作,最后再去掉填充的值。

三、指数加权移动平均平滑

指数加权移动平均平滑算法的原理主要涉及对历史数据进行加权平均,并赋予近期数据更大的权重。

计算方法
EWMA赋予每个数据点的权重随时间呈指数式递减,即越靠近当前时刻的数据点权重越大。权重分配是通过一个平滑系数α来实现的。该系数决定了近期数据相对于历史数据的权重比例。较大的α值意味着当前数据点的权重更大,平滑效果更灵敏于近期的变化;而较小的α值则使得平滑结果更加平滑,但可能会引入一定的滞后性。

实现代码

import pandas as pd
import matplotlib.pyplot as plt  

# 读取数据
csv_file1 = "outputdata.csv"
'''
csv文件中的数据格式为:
x,y
6.15,1.87
6.13,1.54
6.05,1.46
...

'''
df1 = pd.read_csv(csv_file1)

# 指数加权移动平均平滑
span = 5   #   α = 2 / (span + 1)  span越大,α越小
df1['x_smoothed_ewma'] = df1['x'].ewm(span=span).mean()
df1['y_smoothed_ewma'] = df1['y'].ewm(span=span).mean()

plt.figure(figsize=(10, 6))  # 设置图形大小  
plt.plot(df1['x'], df1['y'], linewidth=4, label='original')  # 绘制原始数据  
plt.plot(df1['x_smoothed_ewma'], df1['y_smoothed_ewma'], linewidth=2, label='ma')  # 绘制滤波数据 

# 添加图例  
plt.legend()  
  
# 设置纵横比相等,这取决于数据的具体范围  
plt.gca().set_aspect('equal', adjustable='box')

# 添加标题和坐标轴标签  
plt.title('Data Smoothing')  
plt.xlabel('X Axis')  
plt.ylabel('Y Axis')  
plt.grid(True)

# 显示图形  
plt.show()

平滑效果:
在这里插入图片描述

四、Savitzky-Golay滤波

Savitzky-Golay滤波器通过滑动窗口内的多项式拟合来实现信号的平滑。

计算方法

  1. 选择滤波窗口:首先选择一个合适大小的滑动窗口,这个窗口在数据上滑动,对窗口内的数据进行处理。
  2. 多项式拟合:在每个滑动窗口内,使用多项式函数对数据进行最小二乘法拟合。多项式的阶数和窗口大小可以根据数据的特性进行调整,以达到最佳的滤波效果。
  3. 计算拟合值:根据拟合得到的多项式函数,计算窗口中心点的估计值,作为滤波后的结果。这个估计值反映了窗口内数据的局部趋势,从而实现了数据的平滑。

实现代码

import pandas as pd
from scipy.signal import savgol_filter
import matplotlib.pyplot as plt  

# 读取数据
csv_file1 = "outputdata.csv"
'''
csv文件中的数据格式为:
x,y
6.15,1.87
6.13,1.54
6.05,1.46
...

'''
df1 = pd.read_csv(csv_file1)

window_size = 5  # 窗口大小
# Savitzky-Golay 滤波
df1['x_smoothed_sg'] = savgol_filter(df1['x'], window_length=window_size, polyorder=2)
df1['y_smoothed_sg'] = savgol_filter(df1['y'], window_length=window_size, polyorder=2)

plt.figure(figsize=(10, 6))  # 设置图形大小  
plt.plot(df1['x'], df1['y'], linewidth=4, label='original')  # 绘制原始数据  
plt.plot(df1['x_smoothed_sg'], df1['y_smoothed_sg'], linewidth=2, label='ma')  # 绘制滤波数据 

# 添加图例  
plt.legend()  
  
# 设置纵横比相等,这取决于数据的具体范围  
plt.gca().set_aspect('equal', adjustable='box')

# 添加标题和坐标轴标签  
plt.title('Data Smoothing')  
plt.xlabel('X Axis')  
plt.ylabel('Y Axis')  
plt.grid(True)

# 显示图形  
plt.show()

平滑效果:
在这里插入图片描述

五、高斯滤波

高斯滤波计算滤波窗口内各点相对于中心点的权重。权重分配的原则是:距离中心点越近的点权重越大,距离越远的点权重越小。这样,通过对窗口内各点进行加权平均,可以实现数据的平滑处理。

实现代码

import pandas as pd
from scipy.ndimage import gaussian_filter1d
import matplotlib.pyplot as plt  

# 读取数据
csv_file1 = "outputdata.csv"
'''
csv文件中的数据格式为:
x,y
6.15,1.87
6.13,1.54
6.05,1.46
...

'''
df1 = pd.read_csv(csv_file1)

# 高斯滤波
sigma = 1  # 高斯核的标准差
df1['x_smoothed_gaussian'] = gaussian_filter1d(df1['x'], sigma)
df1['y_smoothed_gaussian'] = gaussian_filter1d(df1['y'], sigma)

plt.figure(figsize=(10, 6))  # 设置图形大小  
plt.plot(df1['x'], df1['y'], linewidth=4, label='original')  # 绘制原始数据  
plt.plot(df1['x_smoothed_gaussian'], df1['y_smoothed_gaussian'], linewidth=2, label='ma')  # 绘制滤波数据 

# 添加图例  
plt.legend()  
  
# 设置纵横比相等,这取决于数据的具体范围  
plt.gca().set_aspect('equal', adjustable='box')

# 添加标题和坐标轴标签  
plt.title('Data Smoothing')  
plt.xlabel('X Axis')  
plt.ylabel('Y Axis')  
plt.grid(True)

# 显示图形  
plt.show()

平滑效果:
在这里插入图片描述

六、卡尔曼滤波

卡尔曼滤波是一种动态系统状态估计的算法,它能够从带有噪声的观测数据中估计出系统的最优状态。卡尔曼滤波基于系统的状态空间模型,通过不断地预测和更新状态来实现对系统状态的跟踪和估计。

下面介绍了卡尔曼算法的一个简单应用,其中:
dim_x:状态向量的维度。在这里设置为 2,表示我们定义的系统状态有两个变量。
dim_z:观测向量的维度。这里设置为 1,因为我们每次观测到的是一个值(即 y 的值)。
F 是状态转移矩阵,描述了系统从一个时刻到下一个时刻的状态如何变化。表示位置在下一个时刻等于当前位置加上速度,速度在下一个时刻保持不变。
H 是观测矩阵,将状态向量与观测值联系起来。这里表示观测值只与状态向量中的位置有关,与速度无关。
x 是初始状态向量,这里将初始位置和速度都设置为 0。
P 是初始协方差矩阵,它表示对初始状态向量的不确定性估计。
np.eye(2) 创建一个 2x2 的单位矩阵,然后乘以 10 表示对位置和速度的初始不确定性都较大。
predict() 方法根据状态转移矩阵 F 预测下一个状态,并更新协方差矩阵 P
update(measurement) 方法根据观测值 measurement 来更新状态向量 x 和协方差矩阵 P

实现代码

import pandas as pd
import numpy as np
from filterpy.kalman import KalmanFilter
import matplotlib.pyplot as plt  

# 读取数据
csv_file1 = "outputdata.csv"
'''
csv文件中的数据格式为:
x,y
6.15,1.87
6.13,1.54
6.05,1.46
...

'''
df1 = pd.read_csv(csv_file1)

# 创建 x 的卡尔曼滤波器
kf_x = KalmanFilter(dim_x=2, dim_z=1)
kf_x.F = np.array([[1, 1], [0, 1]])
kf_x.H = np.array([[1, 0]])
kf_x.x = np.array([[0], [0]])
kf_x.P = np.eye(2) * 10

# 创建 y 的卡尔曼滤波器
kf_y = KalmanFilter(dim_x=2, dim_z=1)
kf_y.F = np.array([[1, 1], [0, 1]])
kf_y.H = np.array([[1, 0]])
kf_y.x = np.array([[0], [0]])
kf_y.P = np.eye(2) * 10

# 进行卡尔曼滤波
for i in range(len(df1)):
    measurement_x = df1['x'].iloc[i]
    kf_x.predict()
    kf_x.update(measurement_x)
    df1.loc[i, 'x_smoothed_kalman'] = kf_x.x[0, 0]

    measurement_y = df1['y'].iloc[i]
    kf_y.predict()
    kf_y.update(measurement_y)
    df1.loc[i, 'y_smoothed_kalman'] = kf_y.x[0, 0]

plt.figure(figsize=(10, 6))  # 设置图形大小  
plt.plot(df1['x'], df1['y'], linewidth=4, label='original')  # 绘制原始数据  
plt.plot(df1['x_smoothed_kalman'], df1['y_smoothed_kalman'], linewidth=2, label='ma')  # 绘制滤波数据 

# 添加图例  
plt.legend()  
  
# 设置纵横比相等,这取决于数据的具体范围  
plt.gca().set_aspect('equal', adjustable='box')

# 添加标题和坐标轴标签  
plt.title('Data Smoothing')  
plt.xlabel('X Axis')  
plt.ylabel('Y Axis')  
plt.grid(True)

# 显示图形  
plt.show()

平滑效果:
在这里插入图片描述

  • 48
    点赞
  • 27
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

Limiiiing

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值