Hampel滤波器作为一种强大的时间序列异常值处理工具，在数据清洗和分析中发挥着重要作用。 本文将深入研究Hampel滤波器的原理和数学推导，并通过实际代码演示其在异常值处理中的应用
1. Hampel滤波器简介
1.1 什么是Hampel滤波器？
Hampel滤波器是一种基于中值和中值绝对偏差（MAD）的滤波器，旨在识别和去除时间序列数据中的异常值。相对于传统均值和标准差方法，Hampel滤波器对异常值更具鲁棒性
1.2 为什么选择Hampel滤波器？
在真实世界的数据中，异常值往往会对统计分析产生较大影响。Hampel滤波器通过使用中值和MAD，适应异常值的存在，提高异常值检测的准确性
2. Hampel滤波器原理
2.1 中值和MAD
Hampel滤波器的核心在于中值的计算和MAD的求解。中值代表数据的中间值，而MAD度量了数据点与中值之间的离散程度
2.2 数学推导
中值：
对于包含 个数据点的数据集 ，中值计算公式如下：
- 如果 为奇数： 中值
- 如果 为偶数： 中值
中值绝对偏差（MAD）：
MAD是每个数据点 与数据集 的中值之间绝对差的中值：
中值中值
其中， 表示每一个数据点
Hampel滤波器异常值判定标准：
如果数据点 被认定为异常值（即其绝对偏差除以MAD超过阈值），那么将其替换为数据集 的中值；否则，保持数据点的原始值
中值中值阀值
3. 代码实例
3.1 数据生成与可视化
首先，我们生成一个包含异常值的正弦波数据，并通过可视化展示原始数据，感受一下异常值在数据集的呈现
import matplotlib.pyplot as plt
import numpy as np
# 设置中文字体
plt.rcParams['font.sans-serif'] = 'SimHei'
plt.rcParams['axes.unicode_minus'] = False
# 生成包含异常值的正弦波数据
def generate_data_with_outliers():
time = np.linspace(0, 10, 100)
signal = np.sin(time) + np.random.normal(0, 0.1, 100)
# 添加异常值
outliers_indices = [20, 40, 60, 80]
outliers_values = [2.0, -1.9, 2.1, -0.5]
for index, value in zip(outliers_indices, outliers_values):
signal[index] = value
return time, signal
# 生成包含异常值的数据
time, data = generate_data_with_outliers()
# 获取异常值的索引和数值
outliers_indices = [20, 40, 60, 80]
outliers_values = [data[i] for i in outliers_indices]
# 绘制原始数据
plt.figure(figsize=(10, 6))
plt.plot(time, data, label='原始数据')
# 仅将异常值标红
plt.scatter(time[outliers_indices], outliers_values, c='red', marker='o', label='异常值', s=50)
plt.title('带有异常值的原始数据')
plt.legend()
plt.show()
3.2 Hampel滤波器定义
接下来，定义Hampel滤波器函数
def hampel(vals_orig, k=7, t0=3):
"""
使用Hampel滤波器去除时间序列中的异常值。
参数：
- vals_orig: numpy数组，原始时间序列数据
- k: 整数，滤波器窗口的大小（半窗口大小为k/2）
- t0: 浮点数，用于异常值检测的阈值
返回：
- vals_filt: numpy数组，经过滤波的时间序列数据
- outliers_indices: list，异常值的索引列表
"""
# 创建输入数据的副本
vals_filt = np.copy(vals_orig)
outliers_indices = []
# 定义Hampel滤波器函数
n = len(vals_orig)
for i in range(k, n - k):
# 提取窗口
window = vals_orig[i - k:i + k + 1]
# 计算中值和中值绝对偏差（MAD）
median = np.median(window)
mad = np.median(np.abs(window - median))
# 检查当前值是否为异常值
if np.abs(vals_orig[i] - median) > t0 * mad:
# 用中值替换异常值
vals_filt[i] = median
# 记录异常值的索引
outliers_indices.append(i)
return vals_filt, outliers_indices
3.3 Hampel滤波器实现
接下来，我们使用Hampel滤波器对数据进行处理，去除异常值
# 对数据应用Hampel滤波器
filtered_data, outliers_indices = hampel(data)
3.4 滤波效果展示
最后，我们展示Hampel滤波器处理后的数据，并通过可视化对比原始数据，特别突出异常值的识别和去除效果
# 绘制原始数据和滤波后的数据
plt.figure(figsize=(10, 6))
plt.plot(time, data, label='原始数据')
plt.scatter(time[outliers_indices], data[outliers_indices], c='red', marker='o', label='异常值', s=50)
plt.plot(time, filtered_data, label='滤波后数据', linestyle='--', color='red')
plt.fill_between(time, filtered_data, data, color='red', alpha=0.2, label='异常值区域')
plt.title('使用Hampel滤波器去除异常值')
plt.grid(True)
plt.legend()
plt.show()
# 打印异常值的索引
print("异常值的索引：", outliers_indices)
结论
Hampel滤波器通过使用中值和MAD，对于异常值的鲁棒性较强。在实际应用中，它在清理包含异常值的时间序列数据方面表现出色，是数据分析中一种强大的工具
如果你对类似于这样的文章感兴趣。
欢迎关注、点赞、转发~