Ctrl+K

规则网格上的多元数据插值 (RegularGridInterpolator)#

假设您在规则网格上拥有 N 维数据，并且想要对其进行插值。在这种情况下，RegularGridInterpolator 很有用。支持几种插值策略：最近邻、线性以及奇数阶的张量积样条。

严格来说，此类可以有效地处理在直角网格上给出的数据：具有可能不均匀点间距的超立方晶格。每个维度的点数对于不同的维度可以不同。

以下示例演示了它的用法，并比较了使用每种方法的插值结果。

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

假设我们要插值这个二维函数。

>>> def F(u, v):
...     return u * np.cos(u * v) + v * np.sin(u * v)

假设我们只知道规则网格上的一些数据。

>>> fit_points = [np.linspace(0, 3, 8), np.linspace(0, 3, 11)]
>>> values = F(*np.meshgrid(*fit_points, indexing='ij'))

创建用于评估的测试点和真实值。

>>> ut, vt = np.meshgrid(np.linspace(0, 3, 80), np.linspace(0, 3, 80), indexing='ij')
>>> true_values = F(ut, vt)
>>> test_points = np.array([ut.ravel(), vt.ravel()]).T

我们可以创建插值器并使用每种方法插值测试点。

>>> interp = RegularGridInterpolator(fit_points, values)
>>> fig, axes = plt.subplots(2, 3, figsize=(10, 6))
>>> axes = axes.ravel()
>>> fig_index = 0
>>> for method in ['linear', 'nearest', 'slinear', 'cubic', 'quintic']:
...     im = interp(test_points, method=method).reshape(80, 80)
...     axes[fig_index].imshow(im)
...     axes[fig_index].set_title(method)
...     axes[fig_index].axis("off")
...     fig_index += 1
>>> axes[fig_index].imshow(true_values)
>>> axes[fig_index].set_title("True values")
>>> fig.tight_layout()
>>> fig.show()

正如预期的那样，较高阶的样条插值最接近真实值，但与 linear 或 nearest 相比，计算成本更高。slinear 插值也与 linear 插值匹配。

../../_images/ND_regular_grid-1.png

如果您的数据使得样条方法产生振铃，您可以考虑使用 method="pchip"，它使用 PCHIP 插值器的张量积，每个维度一个 PchipInterpolator。

如果您喜欢函数式接口而不是显式创建类实例，interpn 便利函数提供了等效的功能。

具体来说，以下两种形式给出相同的结果

>>> from scipy.interpolate import interpn
>>> rgi = RegularGridInterpolator(fit_points, values)
>>> result_rgi = rgi(test_points)

和

>>> result_interpn = interpn(fit_points, values, test_points)
>>> np.allclose(result_rgi, result_interpn, atol=1e-15)
True

对于限制在 N 维空间的 (N-1) 维子空间内的数据，即当其中一个网格轴的长度为 1 时，沿此轴的推断由 fill_value 关键字参数控制。

>>> x = np.array([0, 5, 10])
>>> y = np.array([0])
>>> data = np.array([[0], [5], [10]])
>>> rgi = RegularGridInterpolator((x, y), data,
...                               bounds_error=False, fill_value=None)
>>> rgi([(2, 0), (2, 1), (2, -1)])   # extrapolates the value on the axis
array([2., 2., 2.])
>>> rgi.fill_value = -101
>>> rgi([(2, 0), (2, 1), (2, -1)])
array([2., -101., -101.])

注意

如果输入数据使得输入维度的单位不相称且相差多个数量级，则插值器可能存在数值伪影。考虑在插值之前重新缩放数据。

均匀间隔的数据#

如果您正在处理具有整数坐标的笛卡尔网格上的数据，例如重采样图像数据，则这些例程可能不是最佳选择。考虑改用 scipy.ndimage.map_coordinates。

对于等间距网格上的浮点数据，可以轻松地将 map_coordinates 包装成类似 RegularGridInterpolator 的形式。以下是源自 Johanness Buchner 的 ‘regulargrid’ 包 的一个最简示例。

class CartesianGridInterpolator:
    def __init__(self, points, values, method='linear'):
        self.limits = np.array([[min(x), max(x)] for x in points])
        self.values = np.asarray(values, dtype=float)
        self.order = {'linear': 1, 'cubic': 3, 'quintic': 5}[method]

    def __call__(self, xi):
        """
        `xi` here is an array-like (an array or a list) of points.

        Each "point" is an ndim-dimensional array_like, representing
        the coordinates of a point in ndim-dimensional space.
        """
        # transpose the xi array into the ``map_coordinates`` convention
        # which takes coordinates of a point along columns of a 2D array.
        xi = np.asarray(xi).T

        # convert from data coordinates to pixel coordinates
        ns = self.values.shape
        coords = [(n-1)*(val - lo) / (hi - lo)
                  for val, n, (lo, hi) in zip(xi, ns, self.limits)]

        # interpolate
        return map_coordinates(self.values, coords,
                               order=self.order,
                               cval=np.nan)  # fill_value

此包装器可以用作 RegularGridInterpolator 的（几乎）直接替代品。

>>> x, y = np.arange(5), np.arange(6)
>>> xx, yy = np.meshgrid(x, y, indexing='ij')
>>> values = xx**3 + yy**3
>>> rgi = RegularGridInterpolator((x, y), values, method='linear')
>>> rgi([[1.5, 1.5], [3.5, 2.6]])
array([ 9. , 64.9])
>>> cgi = CartesianGridInterpolator((x, y), values, method='linear')
>>> cgi([[1.5, 1.5], [3.5, 2.6]])
array([ 9. , 64.9])

请注意，上面的示例使用 map_coordinates 边界条件。因此，cubic 和 quintic 插值的结果可能与 RegularGridInterpolator 的结果不同。有关边界条件和其他附加参数的更多详细信息，请参阅 scipy.ndimage.map_coordinates 文档。最后，我们注意到此简化示例假设输入数据是按升序给出的。

在本页中