从传感器噪声到清晰趋势：手把手教你用Python重现经典信号预处理案例（含代码避坑）

从传感器噪声到清晰趋势：Python信号预处理实战指南

传感器数据如同未经雕琢的玉石，原始信号中总是混杂着各种噪声。温度传感器会受环境干扰，电压测量躲不开工频干扰，音频信号更是噪声的重灾区。本文将用Python带领你穿越噪声迷雾，还原信号本真。

1. 温度数据平滑：从小时波动到日趋势

波士顿洛根机场的温度数据集是时间序列分析的经典案例。原始数据包含2011年1月每小时温度读数，呈现出明显的昼夜波动。我们的目标是消除这些高频波动，提取出更有意义的日趋势。

1.1 移动平均滤波器的Python实现

MATLAB中的filter函数在Python中可以用numpy.convolve替代，但需要注意边界处理方式的不同：

import numpy as np import pandas as pd # 加载温度数据（模拟波士顿机场数据） tempC = pd.read_csv('boston_temp.csv')['Temp'].values hours_per_day = 24 coeff_24h = np.ones(hours_per_day) / hours_per_day # 使用卷积实现移动平均 avg_24h_tempC = np.convolve(tempC, coeff_24h, mode='same')

注意：mode='same'确保输出长度与输入相同，但会导致边界效应，这与MATLAB的默认行为不同。

1.2 延迟补偿与可视化

移动平均滤波器会引入(N-1)/2的延迟，这在Python中需要手动补偿：

import matplotlib.pyplot as plt days = np.arange(len(tempC)) / 24 f_delay = (len(coeff_24h) - 1) / 2 plt.figure(figsize=(10, 6)) plt.plot(days, tempC, label='Hourly Temp') plt.plot(days - f_delay/24, avg_24h_tempC, label='24 Hour Average') plt.ylabel('Temp (℃)') plt.xlabel('Time elapsed from Jan 1, 2011 (days)') plt.legend() plt.show()

1.3 昼夜温差分析

通过减去平滑后的趋势，我们可以分析昼夜温差模式：

delta_tempC = tempC - avg_24h_tempC delta_tempC = delta_tempC.reshape(31, 24) # 31天×24小时 mean_daily_diff = np.mean(delta_tempC, axis=0) hours = np.arange(24) plt.figure(figsize=(8, 5)) plt.plot(hours, mean_daily_diff) plt.title('Mean temperature differential from 24h average') plt.xlabel('Hour of day') plt.ylabel('Temperature difference (℃)') plt.grid(True)

2. 工频噪声消除：电力线干扰的克星

60Hz的电力线干扰是传感器测量中的常见问题。下面我们处理一个受60Hz噪声污染的电压信号。

2.1 采样率与滤波器设计的关系

from scipy import signal # 模拟受60Hz干扰的电压信号 fs = 1000 # 采样率1kHz t = np.arange(0, 1, 1/fs) f_noise = 60 # 噪声频率 open_loop_voltage = 5 + 0.1 * np.sin(2*np.pi*f_noise*t) + 0.05*np.random.randn(len(t)) # 设计17点移动平均滤波器 window_size = 17 b = np.ones(window_size)/window_size filtered_voltage = signal.lfilter(b, 1, open_loop_voltage)

2.2 重采样技术提升滤波效果

# 重采样到1020Hz以匹配60Hz周期 fs_resamp = 1020 num_samples = int(len(open_loop_voltage) * fs_resamp / fs) resampled_voltage = signal.resample(open_loop_voltage, num_samples) # 应用移动平均滤波器 filtered_resampled = signal.lfilter(b, 1, resampled_voltage) # 绘制结果对比 plt.figure(figsize=(12, 6)) plt.subplot(211) plt.plot(t, open_loop_voltage, label='Original') plt.plot(t, filtered_voltage, label='Filtered at 1kHz') plt.legend() t_resamp = np.arange(len(resampled_voltage)) / fs_resamp plt.subplot(212) plt.plot(t_resamp, resampled_voltage, label='Resampled') plt.plot(t_resamp, filtered_resampled, label='Filtered at 1020Hz') plt.legend()

3. 高级滤波技术：超越简单平均

3.1 Savitzky-Golay滤波器实战

Python中scipy.signal.savgol_filter函数比MATLAB的sglolayfilt需要更明确的参数：

# 三次多项式，窗口长度7 cubic_ma = signal.savgol_filter(tempC, window_length=7, polyorder=3) # 五次多项式，窗口长度9 quintic_ma = signal.savgol_filter(tempC, window_length=9, polyorder=5) plt.figure(figsize=(10, 6)) plt.plot(days, tempC, label='Hourly') plt.plot(days, cubic_ma, label='Cubic (win=7)') plt.plot(days, quintic_ma, label='Quintic (win=9)') plt.legend()

关键参数选择原则：

• 窗口长度应大于多项式阶数
• 奇数窗口确保对称性
• 较大窗口提供更强平滑但可能丢失细节

3.2 中值滤波与Hampel滤波替代方案

Python生态中没有直接对应的Hampel滤波器，但可以用scipy.stats和pandas组合实现：

from scipy import stats def hampel_filter(x, window_size=5, n_sigmas=3): k = 1.4826 # 高斯分布的比例因子 rolling_median = pd.Series(x).rolling(window=window_size, center=True).median() residuals = np.abs(x - rolling_median) mad = k * pd.Series(residuals).rolling(window=window_size, center=True).median() threshold = n_sigmas * mad outliers = residuals > threshold x_filtered = np.where(outliers, rolling_median, x) return x_filtered # 测试含离群值的信号 y = np.sin(np.linspace(0, 10, 500)) + 0.1*np.random.randn(500) y[::50] += 2 # 添加离群值 y_filtered = hampel_filter(y, window_size=13)

4. 实战避坑指南：Python与MATLAB的关键差异

4.1 边界处理行为对比

MATLAB的smoothdata函数默认使用'zeropad'处理边界，而Python的scipy.signal函数通常提供多种选项：

4.2 滤波器延迟补偿技巧

Python中实现零相位滤波的正确方式：

# 使用filtfilt实现零相位滤波 zero_phase_filtered = signal.filtfilt(b, 1, tempC) # 对比常规滤波 plt.figure(figsize=(10, 6)) plt.plot(days, tempC, label='Original') plt.plot(days, avg_24h_tempC, label='Causal filter') plt.plot(days, zero_phase_filtered, label='Zero-phase filter') plt.legend()

4.3 性能优化策略

处理长信号时的内存优化技巧：

from numba import jit @jit(nopython=True) def moving_average_numba(x, window): result = np.empty_like(x) half = window // 2 for i in range(len(x)): start = max(0, i - half) end = min(len(x), i + half + 1) result[i] = np.mean(x[start:end]) return result # 测试100万点数据 big_data = np.random.randn(1_000_000) %timeit moving_average_numba(big_data, 21) # 比纯Python快10-100倍

信号预处理既是科学也是艺术。在实际项目中，我经常发现没有"最好"的滤波器，只有"最合适"的解决方案。温度数据可能对Savitzky-Golay响应良好，而ECG信号可能更需要小波变换。关键是要理解每种方法的优缺点，并通过可视化不断验证结果是否符合预期。