在二维数组中使用 interp1d 的更快方法

Faster way of using interp1d in 2d array

结果正确。但是在我的实际问题中,数据太大了,所以我想直接应用插值而不使用 循环。任何想法将不胜感激。

import numpy as np
from scipy.interpolate import interp1d


data = np.array([[99,0,3,4,5],
               [6,7,0,9,10],
               [11,22,0,14,15]],dtype=np.float32)

data[data==0] = np.nan

def gap_fill(y):
    not_nan = ~np.isnan(y)
    x = np.arange(len(y))
    interp = interp1d(x[not_nan], y[not_nan], kind='linear')
    ynew = interp(x)
    return ynew

results = []
for d in data:
   gapfilled = gap_fill(d)
   results.append(gapfilled)
print results

[array([ 99.,  51.,   3.,   4.,   5.]), array([  6.,   7.,   8.,   9.,  10.]), array([ 11.,  22.,  18.,  14.,  15.])]

我一时兴起想到的是:

In [8]: gap_fill(data.flatten()).reshape(data.shape)
Out[8]: 
array([[ 99.,  51.,   3.,   4.,   5.],
       [  6.,   7.,   8.,   9.,  10.],
       [ 11.,  22.,  18.,  14.,  15.]])

这适用于您的示例,因为所有 nan 都在行内部。但是,对于行末尾的元素,这会将外推转换为跨行插值,这可能是您不想要的。

严格来说,线性插值是在两点 (1-a)*x1+a*x2 之间找到值,其中 0<=a<=1。如果 a 超出该范围,则为线性外推。

interp1 中的默认操作是在外推情况下引发错误。由于您的迭代 gap_fill 运行,因此您不能有任何外推案例。在这种情况下,我的展平解决方案应该可以正常工作。

看起来interp1d 没有使用任何 C 代码进行线性插值。同时查看其文档,您可能会通过添加 copy=False, assume_sorted=True.

来提高速度

它的核心作用是:

    slope = (y_hi - y_lo) / (x_hi - x_lo)[:, None]
    y_new = slope*(x_new - x_lo)[:, None] + y_lo