两个二维数组的线性插值

Question

在上一个问题（在二维数组上使用 numpy.interp的最快方法）中，有人要求以最快的方式实现以下内容：

np.array([np.interp(X[i], x, Y[i]) for i in range(len(X))])

假设X和Y是具有多行的矩阵，因此 for 循环代价高昂。在这种情况下有一个很好的解决方案可以避免 for 循环（请参阅上面的链接答案）。

---

我面临着一个非常相似的问题，但我不清楚在这种情况下是否可以避免 for 循环：

np.array([np.interp(x, X[i], Y[i]) for i in range(len(X))])

换句话说，我想使用线性插值对存储在两个矩阵X和的行中的大量信号进行上采样Y。我希望在 numpy 或 scipy (scipy.interpolate.interp1d) 中找到一个通过广播语义支持此操作的函数，但到目前为止我似乎找不到一个。

其他要点：

• 如果有帮助，行X[i]和x已在我的应用程序中预先排序。此外，就我而言，len(x)它比len(X[i]).

• 该函数scipy.signal.resample几乎可以满足我的要求，但它不使用线性插值...

Answer 1

这是一种直接实现线性插值的矢量化方法。首先，对于每个 x 值和每个 i，j 计算权重，w表示间隔 (X[i, j], X[i, j+1]) 位于 x 左侧的程度。

• 如果整个区间都在 x 的左侧，则该区间的权重为 1。
• 如果没有子区间位于左侧，则权重为 0
• 否则，权重是 0 到 1 之间的数字，表示该区间在 x 左侧的比例。

然后 PL 插值的值计算为 Y[i, 0] + 差值之和 dY[i, j] 乘以相应的权重。逻辑是遵循插值在不同间隔之间的变化量。差异dY = np.diff(Y, axis=1)显示了它在整个时间间隔内的变化程度。乘以相应变化的重量按比例分配。

设置，带有一些小数据数组

import numpy as np
X = np.array([[0, 2, 5, 6, 9], [1, 3, 4, 7, 8]])
Y = np.array([[3, 5, 2, 4, 1], [8, 6, 9, 5, 4]])
x = np.linspace(1, 8, 20)

计算

dX = np.diff(X, axis=1)
dY = np.diff(Y, axis=1)
w = np.clip((x - X[:, :-1, None])/dX[:, :, None], 0, 1)
y = Y[:, [0]] + np.sum(w*dY[:, :, None], axis=1)

示范

这只是为了表明插值是正确的。蓝点：原始数据，红点是计算所得。

import matplotlib.pyplot as plt
plt.plot(x, y[0], 'ro')
plt.plot(X[0], Y[0], 'bo')
plt.plot(x, y[1], 'rd')
plt.plot(X[1], Y[1], 'bd')
plt.show()