如何以正确的方式平滑曲线？

Question

假设我们有一个数据集

import numpy as np
x = np.linspace(0,2*np.pi,100)
y = np.sin(x) + np.random.random(100) * 0.2

因此，我们有20%的数据集变化。我的第一个想法是使用scipy的单变量pline函数，但问题是这并不能很好地考虑小噪声。如果你考虑频率，背景比信号小得多，所以一个想法可能是仅限截止的样条曲线，但这将涉及前后傅里叶变换，这可能会导致不良行为。 另一种方法是移动平均线，但这也需要正确选择延迟时间。

有什么解决这个问题的提示/书籍或链接吗？

Answer 1

如果您对周期性信号的“平滑”版本感兴趣（如您的示例），那么FFT是正确的方法。取傅里叶变换并减去低频：

import numpy as np
import scipy.fftpack

N = 100
x = np.linspace(0,2*np.pi,N)
y = np.sin(x) + np.random.random(N) * 0.2

w = scipy.fftpack.rfft(y)
f = scipy.fftpack.rfftfreq(N, x[1]-x[0])
spectrum = w**2

cutoff_idx = spectrum < (spectrum.max()/5)
w2 = w.copy()
w2[cutoff_idx] = 0

y2 = scipy.fftpack.irfft(w2)

即使你的信号不是完全周期性的，这也能很好地去除白噪声。有许多类型的过滤器可供使用（高通、低通等），合适的过滤器取决于您要找的是什么。

Answer 2

我喜欢Savitzky-Golay filter。它使用最小二乘法将数据的一个小窗口回归到多项式上，然后使用多项式估计窗口中心的点。最后，窗口向前移动一个数据点，过程重复。这种情况一直持续到每一个点相对于它的邻域进行了最佳调整。即使对来自非周期性和非线性源的噪声样本，它也能很好地工作。

这是一张thorough cookbook example。请参阅下面的代码，了解它的易用性。注意：我省略了定义savitzky_golay()函数的代码，因为您可以从上面链接的菜谱示例中复制/粘贴它。

import numpy as np
import matplotlib.pyplot as plt

x = np.linspace(0,2*np.pi,100)
y = np.sin(x) + np.random.random(100) * 0.2
yhat = savitzky_golay(y, 51, 3) # window size 51, polynomial order 3

plt.plot(x,y)
plt.plot(x,yhat, color='red')
plt.show()

更新：我注意到我链接到的食谱示例已被删除。幸运的是，正如@dodohjk所指出的，Savitzky Golay过滤器已经被并入into the SciPy library。 要使用SciPy源代码改编上述代码，请键入：

from scipy.signal import savgol_filter
yhat = savgol_filter(y, 51, 3) # window size 51, polynomial order 3

Answer 3

基于移动平均值框（通过卷积）的一种快速而肮脏的平滑数据方法：

x = np.linspace(0,2*np.pi,100)
y = np.sin(x) + np.random.random(100) * 0.8

def smooth(y, box_pts):
    box = np.ones(box_pts)/box_pts
    y_smooth = np.convolve(y, box, mode='same')
    return y_smooth

plot(x, y,'o')
plot(x, smooth(y,3), 'r-', lw=2)
plot(x, smooth(y,19), 'g-', lw=2)