这是我每0.1秒就有一辆车的位置数据的(x,y)曲线图.总分在500分左右.

original data

我读到了其他关于使用SciPy(herehere)进行内插的解决方案,但似乎在默认情况下,SciPy以均匀的间隔进行内插.以下是我当前的代码:

def reduce_dataset(x_list, y_list, num_interpolation_points):
    points = np.array([x_list, y_list]).T 
    distance = np.cumsum( np.sqrt(np.sum( np.diff(points, axis=0)**2, axis=1 )) )
    distance = np.insert(distance, 0, 0)/distance[-1]
    interpolator =  interp1d(distance, points, kind='quadratic', axis=0)
    results = interpolator(np.linspace(0, 1, num_interpolation_points)).T.tolist()
    new_xs = results[0]
    new_ys = results[1]
    return new_xs, new_ys



xs, ys = reduce_dataset(xs,ys, 50)
colors = cm.rainbow(np.linspace(0, 1, len(ys)))
i = 0
for y, c in zip(ys, colors):
    plt.scatter(xs[i], y, color=c)
    i += 1

它产生以下输出:

interpolator

这是不错的,但我想设置插值器,try 在最难进行线性插补的地方放置更多的点,而在可以使用插值线轻松重建的区域放置较少的点.

请注意,在第二张图片中,最后一个点似乎突然从前一个点"跳"了出来.中间的部分似乎有点多余,因为其中许多点都落在一条完全直线上.对于要使用线性内插法尽可能准确地重建的东西,这不是使用50个点的最有效使用.

我是手工制作的,但我正在寻找类似这样的东西,其中的算法足够智能,可以在数据非线性变化的地方非常密集地放置点:

enter image description here

通过这种方式,可以以更高的准确度对数据进行插值.该图中点之间的大间隙可以用简单的线非常精确地插值,而密集的聚类需要更频繁的采样. 我已经阅读到interpolator docs on SciPy,但似乎找不到任何发电机或设置可以做到这一点.

我也try 过使用"线性"和"三次"插值法,但它似乎仍然以均匀的间隔进行采样,而不是对最需要它们的点进行分组.

这是SciPy可以做到的吗,或者我应该为这样的工作使用类似SKLearn ML算法的东西吗?

推荐答案

在我看来,您混淆了由interp1d构造的插补器对象和实际的插补坐标,后者是您想要的最终结果.

似乎在默认情况下,SciPy以均匀的间隔进行内插

interp1d将返回根据您提供的xy坐标构建的插补器对象.它们根本不需要均匀分布.

然后,您向该模块提供xnew个值,这些值定义了该模块将在何处重建信号.这是你必须指定是否要均匀间隔的地方:results = interpolator(np.linspace(0, 1, num_interpolation_points)).T.tolist().注意对np.linspace的调用,它的字面意思是"线性间隔的值".

将其替换为np.logspace(),以获得几何间隔值,或者替换为其他值:

import numpy as np
from scipy.interpolate import interp1d

import matplotlib.pyplot as plt

# Generate fake data
x = np.linspace(1, 3, 1000)
y = (x - 2)**3

# interpolation
interpolator = interp1d(x, y)

# different xnews
N = 20
xnew_linspace = np.linspace(x.min(), x.max(), N)  # linearly spaced
xnew_logspace = np.logspace(np.log10(x.min()), np.log10(x.max()), N)  # log spaced

# spacing based on curvature
gradient = np.gradient(y, x)
second_gradient = np.gradient(gradient, x)
curvature = np.abs(second_gradient) / (1 + gradient**2)**(3 / 2)
idx = np.round(np.linspace(0, len(curvature) - 1, N)).astype(int)
epsilon = 1e-1

a = (0.99 * x.max() - x.min()) / np.sum(1 / (curvature[idx] + epsilon))
xnew_curvature = np.insert(x.min() + np.cumsum(a / (curvature[idx] + epsilon)), 0, x.min())

fig, axarr = plt.subplots(2, 2, layout='constrained', sharex=True, sharey=True)

axarr[0, 0].plot(x, y)
for ax, xnew in zip(axarr.flatten()[1:], [xnew_linspace, xnew_logspace, xnew_curvature]):
    ax.plot(xnew, interpolator(xnew), '.--')

axarr[0, 0].set_title('base signal')
axarr[0, 1].set_title('linearly spaced')
axarr[1, 0].set_title('log spaced')
axarr[1, 1].set_title('curvature based spaced')

plt.savefig('test_interp1d.png', dpi=400)

enter image description here

请注意,我不确定像我所做的那样对曲率进行zoom 是正确的方法.但这给了你大约interp1d个的 idea .

Python相关问答推荐

在行数据为向量的DataFrame上计算逐行更改

Df.Drop_Duplates(),以极点表示?

如何在JAX中训练具有多输出(向量值)损失函数的梯度下降模型?

将新列添加到现有极点数据框中

如何在Polars DataFrame中使用`isin‘?

SqlalChemy-同时过滤父对象和子对象

TypeError:无法实例化抽象模型(Django)

没有目标的验证数据

如何从标准输入中读取超过4096个字节,并将其复制粘贴到Linux上的终端?

如何使用GSILD获得正确的日期

如何用Set_Layout_Engine替换Set_Tell_Layout?

极点过滤多条件数据帧

当鼠标悬停和/或点击数据点时, map 光标会显示数据

从文本库中的屏幕内部切换模式--Python

嵌套类型规范中的TypeVar

有没有办法只将两个NumPy数组的特定元素相加在一起?

极点:在同一分组上执行多个操作

NumPy能取代这些列表理解,让它运行得更快吗?

合并条目顺序不匹配的两个dataFrame

如果A的子串在一个DF的B的串中,则不 Select 行