对于包含 nan 的数组的行，均值为 25 个百分点

Question

我有一个二维数组 x，每一行都有不同数量的 nan 值：

array([[   nan, -0.355, -0.036, ...,    nan,    nan],
       [   nan, -0.341, -0.047, ...,    nan,  0.654],
       [  .016, -1.147, -0.667, ...,    nan,    nan],
       ..., 
       [   nan,  0.294, -0.235, ...,    0.65,   nan]])

鉴于此数组，对于每一行，我想计算前 25 个百分位数内所有值的平均值。我正在执行以下操作：

limit = np.nanpercentile(x, 25, axis=1) # output 1D array
ans = np.nanmean(x * (x < limit[:,None]), axis=1)

但这给出了错误的结果——特别是计数 (np.nansum/np.nanmean) 保持不变，无论我选择什么百分位，因为比较在不正确的地方产生零，并被计为均值的有效值。我不能简单地使用 x[x>limit[:,None]] 因为它给出了一个一维数组，我需要一个二维结果。

我通过以下方式解决了它：

f = x.copy()
f[f > limit[:,None]] = np.nan
ans = np.nanmean(f, axis=1)

有更好的方法吗？

Answer 1

方法 #1： 您可以创建一个无效掩码，即来自原始数组的 NaNs 和来自 f > limit[:,None] 的掩码。然后，使用此掩码执行 np.nanmean 等效方法，只考虑 masking 的有效方法。使用 masks/boolean arrays 的好处在于内存方面，因为它占用的内存比浮动 pt 数组少 8 倍。因此，我们将有一个像这样的实现 -

# Create mask of non-NaNs and thresholded ones
mask = ~np.isnan(x) & (x <= limit[:,None])

# Get the row, col indices. Use the row indices for bin-based summing and
# finally averaging by using those indices to get the group lengths.
r,c = np.where(mask)
out = np.bincount(r,x[mask])/np.bincount(r)

方法 #2： 我们也可以使用 np.add.reduceat，这在这里很有用，因为 bin 已经按照掩码进行了排序。所以，效率会更高一些 -

# Get the valid mask as before
mask = ~np.isnan(x) & (x <= limit[:,None])

# Get valid row count. Use np.add.reduceat to perform grouped summations
# at intervals separated by row indices.
rowc = mask.sum(1)
out = np.add.reduceat(x[mask],np.append(0,rowc[:-1].cumsum()))/rowc

基准测试

函数定义 -

def original_app(x, limit):
    f = x.copy()
    f[f > limit[:,None]] = np.nan
    ans = np.nanmean(f, axis=1) 
    return ans

def proposed1_app(x, limit):
    mask = ~np.isnan(x) & (x <= limit[:,None])
    r,c = np.where(mask)
    out = np.bincount(r,x[mask])/np.bincount(r)
    return out

def proposed2_app(x, limit):
    mask = ~np.isnan(x) & (x <= limit[:,None])
    rowc = mask.sum(1)
    out = np.add.reduceat(x[mask],np.append(0,rowc[:-1].cumsum()))/rowc
    return out

时间和验证 -

In [402]: # Setup inputs
     ...: x = np.random.randn(400,500)
     ...: x.ravel()[np.random.randint(0,x.size,x.size//4)] = np.nan # Half as NaNs
     ...: limit = np.nanpercentile(x, 25, axis=1)
     ...: 

In [403]: np.allclose(original_app(x, limit),proposed1_app(x, limit))
Out[403]: True

In [404]: np.allclose(original_app(x, limit),proposed2_app(x, limit))
Out[404]: True

In [405]: %timeit original_app(x, limit)
100 loops, best of 3: 5 ms per loop

In [406]: %timeit proposed1_app(x, limit)
100 loops, best of 3: 4.02 ms per loop

In [407]: %timeit proposed2_app(x, limit)
100 loops, best of 3: 2.18 ms per loop

对于包含 nan 的数组的行，均值为 25 个百分点

mean till 25 percentile for rows of an array containing nan

numpy

vectorization

nan

multidimensional-array

array-broadcasting