如何在python中平滑曲线

Posted 2023-03-12

【中文标题】如何在python中平滑曲线

【英文标题】：how to smooth a curve in python

【发布时间】：2014-05-24 05:01:46

【问题描述】：

我有一条熵曲线（1d numpy 数组），但这条曲线有很多噪音。 我想通过平滑删除噪音。

这是我的曲线图：

我已尝试使用 Kaiser-Bessel 滤波器制作卷积产品来解决此问题：

gaussian_curve = window_kaiser(windowLength, beta=20)  # kaiser filter
gaussian_curve = gaussian_curve / sum(gaussian_curve)

for i in range(0, windows_number):
     start = (i * step) + 1
     end = (i * step) + windowLength
     convolution[i] = (np.convolve(entropy[start:end + 1], gaussian_curve, mode='valid'))
     entropy[i] = convolution[i][0]

但此代码返回此错误：

File "/usr/lib/python2.7/dist-packages/numpy/core/numeric.py", line 822, in convolve
    raise ValueError('v cannot be empty')
ValueError: v cannot be empty

具有“有效”模式的 numpy.convolve 运算符返回重叠中的中心元素，但在这种情况下，返回一个空元素。

有没有一种简单的方法来应用平滑？

谢谢！

【问题讨论】：

How to smooth a curve in the right way?的可能重复

【参考方案1】：

好的，我解决了。 我使用了另一种方法：Savitzky-Golay filter

代码：

def savitzky_golay(y, window_size, order, deriv=0, rate=1):

    import numpy as np
    from math import factorial

    try:
        window_size = np.abs(np.int(window_size))
        order = np.abs(np.int(order))
    except ValueError, msg:
        raise ValueError("window_size and order have to be of type int")
    if window_size % 2 != 1 or window_size < 1:
        raise TypeError("window_size size must be a positive odd number")
    if window_size < order + 2:
        raise TypeError("window_size is too small for the polynomials order")
    order_range = range(order+1)
    half_window = (window_size -1) // 2
    # precompute coefficients
    b = np.mat([[k**i for i in order_range] for k in range(-half_window, half_window+1)])
    m = np.linalg.pinv(b).A[deriv] * rate**deriv * factorial(deriv)
    # pad the signal at the extremes with
    # values taken from the signal itself
    firstvals = y[0] - np.abs( y[1:half_window+1][::-1] - y[0] )
    lastvals = y[-1] + np.abs(y[-half_window-1:-1][::-1] - y[-1])
    y = np.concatenate((firstvals, y, lastvals))
    return np.convolve( m[::-1], y, mode='valid')

现在，我可以输入：

entropy = np.array(entropy)
entropy = savitzky_golay(entropy, 51, 3) # window size 51, polynomial order 3

结果是这样的：

【讨论】：

仅供参考：scipy 0.14.0（即将发布）在scipy.signal.savgol_filter 中实现了 Savitzky-Golay 过滤器。