使用scipy在第一轴上插值数组

Question

我在尝试对一个大小为 Nx*Ny*Nz 的数据立方体进行简单的线性插值，想要生成一个新的立方体，同时保持另外两个维度不变，也就是说，最终的输出应该是 Nxnew*Ny*Nz。看起来 scipy 里的 RegularGridInterpolator 是个不错的选择，不过对我来说，如何生成输入数据还是不太明白。

from scipy.interpolate import RegularGridInterpolator
import numpy as np

x = np.linspace(1,4,11)
y = np.linspace(4,7,22)
z = np.linspace(7,9,33)
V = np.zeros((11,22,33))
for i in range(11):
    for j in range(22):
        for k in range(33):
            V[i,j,k] = 100*x[i] + 10*y[j] + z[k]
fn = RegularGridInterpolator((x,y,z), V)
pts = np.array([[[[2,6,8],[3,5,7]], [[2,6,8],[3,5,7]]]])
out = fn(pts)
print(out, out.shape)

在这个简单的例子中，我想使用新的点 xnew = np.linspace(2,3,50)，同时保持 y 和 z 不变，这样得到的数组形状就变成了 (50,22,33)。另外，如何将这个方法推广到在一个维度上对 n 维数组进行插值，同时保持其他坐标不变呢？

Answer 1

正如评论中提到的，你可以用一个叫做 np.meshgrid 的函数来替代三层嵌套的循环，这样代码会更容易读懂，也更高效。

x, y, z = np.meshgrid(x, y, z, indexing='ij')
V = 100*x + 10*y + z

关于如何生成输入给你的 fn 对象，注意它的 __call__ 方法 需要的输入形状是 (..., ndim)。在这里，... 是你想要的形状（50, 22, 33），而 ndim 是坐标的数量（对于 x、y 和 z 来说是3）。我们可以用 meshgrid 来生成三个单独的数组来表示这些坐标，但为了形成输入给 fn，我们需要把它们组合成一种方式，使得对应坐标的轴放在最后。实现这个的方式有很多，但最简单的方法是使用 np.stack。

from scipy.interpolate import RegularGridInterpolator
import numpy as np

x0 = np.linspace(1, 4, 11)
y0 = np.linspace(4, 7, 22)
z0 = np.linspace(7, 9, 33)

x, y, z = np.meshgrid(x0, y0, z0, indexing='ij')
V = 100*x + 10*y + z

fn = RegularGridInterpolator((x0, y0, z0), V)

xnew = np.linspace(2, 3, 50)
x, y, z = np.meshgrid(xnew, y0, z0, indexing='ij')

xi = np.stack((x, y, z), axis=-1)
# or
# xi = np.moveaxis(np.asarray([x, y, z]), 0, -1)
# or
# xi = np.concatenate((x[..., np.newaxis, ], y[..., np.newaxis], z[..., np.newaxis]), axis=-1)

out = fn(xi)
print(out.shape)
# (50, 22, 33)

你提到的关于将问题推广到“n维数组”的“n”的意思可能有点模糊。一般来说，“n”代表坐标的数量，也就是你的 V 的维度。如果是这样的话，推广到任意数量坐标的情况其实很简单：只需像处理 x、y 和 z 一样处理这些坐标（例如 t、u、v、w）。