【问题标题】:Utilising Savitzky-Golay Filter in R vs Python在 R 与 Python 中使用 Savitzky-Golay 过滤器
【发布时间】:2021-08-16 15:11:23
【问题描述】:

我目前正在尝试在R 中呈现与在Python 中相同的结果,但我认为我一定误解了 Savitzky-Golay 过滤器。我有以下Python 代码:

import numpy as np
from scipy.signal import savgol_filter
t = np.linspace(0,1,10)
X = np.vstack((np.sin(t),np.cos(t))).T
sfd = savgol_filter(X, window_length=5, polyorder=3, axis=0)
sfd
array([[-4.78900581e-07,  9.99997881e-01],
       [ 1.10884544e-01,  9.93841986e-01],
       [ 2.20394870e-01,  9.75397369e-01],
       [ 3.27190431e-01,  9.44944627e-01],
       [ 4.29950758e-01,  9.02837899e-01],
       [ 5.27408510e-01,  8.49596486e-01],
       [ 6.18361741e-01,  7.85877015e-01],
       [ 7.01688728e-01,  7.12465336e-01],
       [ 7.76378020e-01,  6.30281243e-01],
       [ 8.41469460e-01,  5.40300758e-01]])

据我了解,这可以平滑矩阵并准备好开发导数项。但是,当在 R(Savitzky-Golay 函数的最新更新版本)中使用 pracma 时,我得到:

library(pracma)
t = seq(0, 1,length = 10)
X = t(rbind(sin(t), cos(t)))
savgol(X[, 1], fl = 5)
[1] 1.229175e-16 1.108826e-01 2.203977e-01 3.271947e-01 4.299564e-01 5.274154e-01 6.183698e-01 7.016979e-01 7.763719e-01 8.414710e-01

有谁知道为什么这些数字如此不同,以及如何从PythonR 中产生相同的结果?

提前致谢。

【问题讨论】:

    标签: python r scipy differentiation


    【解决方案1】:

    使用信号包中的sgolayfilt函数:

    library(signal)
    packageVersion("signal")
    ## [1] ‘0.7.7’
    
    apply(X, 2, sgolayfilt, n = 5)
    ##                [,1]      [,2]
    ##  [1,] -4.789006e-07 0.9999979
    ##  [2,]  1.108845e-01 0.9938420
    ##  [3,]  2.203949e-01 0.9753974
    ##  [4,]  3.271904e-01 0.9449446
    ##  [5,]  4.299508e-01 0.9028379
    ##  [6,]  5.274085e-01 0.8495965
    ##  [7,]  6.183617e-01 0.7858770
    ##  [8,]  7.016887e-01 0.7124653
    ##  [9,]  7.763780e-01 0.6302812
    ## [10,]  8.414695e-01 0.5403008
    

    【讨论】:

    • 谢谢,我能够做到这一点,但我担心signal 包似乎没有更新,这就是我想进一步了解pracma 功能的原因。
    • 前几天才更新。确保您使用的是答案中的版本。
    • 谢谢,你知道为什么我在尝试上传最新版本时会出错吗?
    • 尝试不同的镜像。也许它还没有达到你正在使用的那个。它也是新的,因此即使您在 Windows 上,您也可能需要从源代码编译它或等到 CRAN 生成二进制文件。如果一切都失败了,直接从它的 CRAN 页面获取它:cran.r-project.org/package=signal
    • 啊,好吧,我在 Mac 上运行,所以可能还没有更新。
    【解决方案2】:

    SciPy 函数savgol_filter 有几个选项用于处理输入数组的末端;请参阅文档字符串中的mode 参数。

    看起来 R 函数 savgol 的行为对应于 SciPy 的 savgol_filter 中的 mode='constant'。除了第一个值(在这两种情况下实际上都是 0)之外,savgol_filter 的输出与 R 中 savgol 的输出相匹配:

    In [82]: sfd = savgol_filter(X, window_length=5, polyorder=4, axis=0, mode='constant')
    
    In [83]: sfd[:, 0]
    Out[83]: 
    array([1.95316193e-17, 1.10882629e-01, 2.20397743e-01, 3.27194697e-01,
           4.29956364e-01, 5.27415386e-01, 6.18369803e-01, 7.01697876e-01,
           7.76371921e-01, 8.41470985e-01])
    

    【讨论】:

    • 啊,好吧,您知道如何在不更改 SciPy 方法的情况下使用 savgol 获得相同的结果吗? R 文档提供 rdocumentation.org/packages/signal/versions/0.7-6/topics/… 这个函数作为“另见”选项,它似乎给出了正确的结果,但我不确定我理解为什么。
    • 如接受的答案中所述,看起来 R signal 包中的 sgolayfiltscipy.signal.savgol_filter 的默认行为匹配。
    猜你喜欢
    • 2016-08-27
    • 1970-01-01
    • 2022-01-13
    • 1970-01-01
    • 1970-01-01
    • 2017-01-29
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多