python/scipy中的多元样条插值？

Question

在python中，是否有库模块或其他简单的方法来实现多元样条插值？

具体来说，我在一个规则间隔的三维网格上有一组标量数据，我需要在分散在整个域中的少量点处进行插值。对于二维，我一直在使用scipy.interpolate.RectBivariateSpline，实际上我在寻找三维数据的扩展。

我发现的N维插值程序还不够好：我更喜欢样条曲线而不是LinearNDInterpolator的平滑度，而且我有太多的数据点（通常超过一百万个）来工作，例如径向基函数。

如果有人知道一个python库可以做到这一点，或者我可以调用或移植另一种语言的python库，我将非常感激。

Answer 1

dim>；2中的平滑样条插值很难实现，因此没有多少免费可用的库能够做到这一点（事实上，我不知道任何库）。

您可以尝试反向距离加权插值，请参见：Inverse Distance Weighted (IDW) Interpolation with Python。 这将产生相当平滑的结果，并比RBF更适合于更大的数据集。

Answer 2

如果我正确理解了你的问题，你输入的“观察”数据会定期被网格化吗？

如果是的话，^{}正是你想要的。

在第一次扫描时有点难以理解，但实际上，你只需要给它一个坐标序列，你想在像素/体素/n维索引坐标中插值网格的值。

作为二维示例：

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

# Note that the output interpolated coords will be the same dtype as your input
# data.  If we have an array of ints, and we want floating point precision in
# the output interpolated points, we need to cast the array as floats
data = np.arange(40).reshape((8,5)).astype(np.float)

# I'm writing these as row, column pairs for clarity...
coords = np.array([[1.2, 3.5], [6.7, 2.5], [7.9, 3.5], [3.5, 3.5]])
# However, map_coordinates expects the transpose of this
coords = coords.T

# The "mode" kwarg here just controls how the boundaries are treated
# mode='nearest' is _not_ nearest neighbor interpolation, it just uses the
# value of the nearest cell if the point lies outside the grid.  The default is
# to treat the values outside the grid as zero, which can cause some edge
# effects if you're interpolating points near the edge
# The "order" kwarg controls the order of the splines used. The default is 
# cubic splines, order=3
zi = ndimage.map_coordinates(data, coords, order=3, mode='nearest')

row, column = coords
nrows, ncols = data.shape
im = plt.imshow(data, interpolation='nearest', extent=[0, ncols, nrows, 0])
plt.colorbar(im)
plt.scatter(column, row, c=zi, vmin=data.min(), vmax=data.max())
for r, c, z in zip(row, column, zi):
    plt.annotate('%0.3f' % z, (c,r), xytext=(-10,10), textcoords='offset points',
            arrowprops=dict(arrowstyle='->'), ha='right')
plt.show()

要在n维中执行此操作，只需传入适当大小的数组：

import numpy as np
from scipy import ndimage

data = np.arange(3*5*9).reshape((3,5,9)).astype(np.float)
coords = np.array([[1.2, 3.5, 7.8], [0.5, 0.5, 6.8]])
zi = ndimage.map_coordinates(data, coords.T)

就缩放和内存使用而言，map_coordinates将创建数组的筛选副本，如果您使用的是order>；1（即不是线性插值）。如果只想在极少数点处插值，则开销相当大。但是，它不会随着要插值的点数而增加。只要有足够的RAM用于输入数据数组的一个临时副本，就可以了。

如果无法将数据副本存储在内存中，可以a）指定prefilter=False和order=1并使用线性插值，或者b）使用ndimage.spline_filter将原始数据替换为筛选版本，然后使用prefilter=False调用map_坐标。

即使您有足够的ram，如果您需要多次调用map_坐标（例如交互使用等），保持过滤后的数据集也会大大加快速度。