python - 具有 scipy.interpolate : Poor interpolation for low-amplitude, 快速振荡函数的样条表示

Question

我需要(以数字方式)计算函数的一阶和二阶导数，为此我尝试同时使用 splrep 和 UnivariateSpline 来创建样条曲线插值函数的导数。

但是，对于幅度为 10^-1 或更低的函数，样条表示本身似乎存在固有问题并且(快速)振荡。

例如，考虑以下代码以在区间 (0,6*pi) 内创建正弦函数的样条表示(因此该函数仅振荡三次):

import scipy
from scipy import interpolate
import numpy
from numpy import linspace
import math
from math import sin

k = linspace(0, 6.*pi, num=10000) #interval (0,6*pi) in 10'000 steps
y=[]
A = 1.e0 # Amplitude of sine function

for i in range(len(k)):

  y.append(A*sin(k[i]))

tck =interpolate.UnivariateSpline(x, y, w=None, bbox=[None, None], k=5, s=2)
M=tck(k)

下面是 A = 1.e0 和 A = 1.e-2 时 M 的结果

http://i.imgur.com/uEIxq.png振幅 = 1

http://i.imgur.com/zFfK0.png振幅 = 1/100

显然，由样条曲线创建的插值函数是完全不正确的!第二张图甚至没有以正确的频率振荡。

有没有人对这个问题有任何见解？或者知道在 numpy/scipy 中创建样条的另一种方法？

干杯，罗里

Answer 1

<罢工>我猜你的问题是由于别名造成的。

什么是 x在你的例子中？

如果 x您插值的值比原始点的间距小，您将固有地丢失频率信息。这完全独立于任何类型的插值。它是下采样所固有的。

不要在意上面关于别名的部分。它不适用于这种情况(尽管我仍然不知道您的示例中的 x 是什么...

我刚刚意识到当您使用非零平滑因子 ( s) 时，您是在 原始输入点上评估您的点。

根据定义，平滑不会完全适合数据。尝试输入 s=0相反。

举个简单的例子:

import matplotlib.pyplot as plt
import numpy as np
from scipy import interpolate

x = np.linspace(0, 6.*np.pi, num=100) #interval (0,6*pi) in 10'000 steps
A = 1.e-4 # Amplitude of sine function
y = A*np.sin(x)

fig, axes = plt.subplots(nrows=2)
for ax, s, title in zip(axes, [2, 0], ['With', 'Without']):
    yinterp = interpolate.UnivariateSpline(x, y, s=s)(x)
    ax.plot(x, yinterp, label='Interpolated')
    ax.plot(x, y, 'bo',label='Original')
    ax.legend()
    ax.set_title(title + ' Smoothing')

plt.show()

您只能清楚地看到低振幅平滑效果的原因是平滑因子的定义方式。请参阅 scipy.interpolate.UnivariateSpline 的文档更多细节。

即使振幅更高，如果使用平滑，插值数据也不会与原始数据匹配。

例如，如果我们只是将振幅 ( A ) 更改为 1.0在上面的代码示例中，我们仍然会看到平滑的效果...