让我们为 python 的 numpy 制作一个 N 维像素分箱/桶化的参考实现



我经常想将像素桶/像素桶作为一个numpy数组,这意味着用单个像素替换N连续像素组,该像素是N替换像素的总和。例如,从以下值开始:

x = np.array([1, 3, 7, 3, 2, 9])

如果bucket大小为2,则转换为:

bucket(x, bucket_size=2) 
= [1+3, 7+3, 2+9]
= [4, 10, 11]

据我所知,没有numpy函数专门做这件事(如果我错了,请纠正我!),所以我经常自己滚动。对于1d numpy数组,这还不错:

import numpy as np
def bucket(x, bucket_size):
    return x.reshape(x.size // bucket_size, bucket_size).sum(axis=1)
bucket_me = np.array([3, 4, 5, 5, 1, 3, 2, 3])
print(bucket(bucket_me, bucket_size=2)) #[ 7 10  4  5]

然而,我很容易对多维情况感到困惑,最终我一次又一次地推出了我自己的错误的、半途而废的解决方案来解决这个"简单"的问题。如果我们能建立一个很好的N维参考实现,我会很高兴的。

  • 优选地,函数调用将允许沿着不同轴的不同仓大小(可能类似于bucket(x, bucket_size=(2, 2, 3))

  • 优选地,该解决方案将是合理有效的(在numpy中整形和求和相当快)

  • 当数组不能很好地划分为整数个桶时,处理边缘效果的加分。

  • 允许用户选择初始垃圾箱边缘偏移的奖励点数。

正如Divakar所建议的,以下是我在一个二维示例中想要的行为:

x = np.array([[1, 2, 3, 4],
              [2, 3, 7, 9],
              [8, 9, 1, 0],
              [0, 0, 3, 4]])
bucket(x, bucket_size=(2, 2))
= [[1 + 2 + 2 + 3, 3 + 4 + 7 + 9],
   [8 + 9 + 0 + 0, 1 + 0 + 3 + 4]]
= [[8, 23],
   [17, 8]]

希望我的算术做对了;)

我认为您可以使用skimage的view_as_blocks来完成大部分繁琐的工作。这个函数是用as_strided实现的,所以它非常有效(它只是改变步幅信息来重塑数组)。因为它是用Python/NumPy编写的,所以如果没有安装skimage,您总是可以复制代码。

应用该函数后,您只需要对整形数组的N个拖尾轴求和(其中N是bucket_size元组的长度)。这里有一个新的bucket()函数:

from skimage.util import view_as_blocks
def bucket(x, bucket_size):
    blocks = view_as_blocks(x, bucket_size)
    tup = tuple(range(-len(bucket_size), 0))
    return blocks.sum(axis=tup)

例如:

>>> x = np.array([1, 3, 7, 3, 2, 9])
>>> bucket(x, bucket_size=(2,))
array([ 4, 10, 11])
>>> x = np.array([[1, 2, 3, 4],
                  [2, 3, 7, 9],
                  [8, 9, 1, 0],
                  [0, 0, 3, 4]])
>>> bucket(x, bucket_size=(2, 2))
array([[ 8, 23],
       [17,  8]])
>>> y = np.arange(6*6*6).reshape(6,6,6)
>>> bucket(y, bucket_size=(2, 2, 3))
array([[[ 264,  300],
        [ 408,  444],
        [ 552,  588]],
       [[1128, 1164],
        [1272, 1308],
        [1416, 1452]],
       [[1992, 2028],
        [2136, 2172],
        [2280, 2316]]])

从as_strided:本地

x = array([[1, 2, 3, 4],
           [2, 3, 7, 9],
           [8, 9, 1, 0],
           [0, 0, 3, 4]])
from numpy.lib.stride_tricks import as_strided     
def bucket(x,bucket_size):
      x=np.ascontiguousarray(x)
      oldshape=array(x.shape)
      newshape=concatenate((oldshape//bucket_size,bucket_size))
      oldstrides=array(x.strides)
      newstrides=concatenate((oldstrides*bucket_size,oldstrides))
      axis=tuple(range(x.ndim,2*x.ndim))
      return as_strided (x,newshape,newstrides).sum(axis)

如果一个维度没有均匀地划分为x的相应维度,则剩余元素将丢失。

验证:

In [9]: bucket(x,(2,2))
Out[9]: 
array([[ 8, 23],
       [17,  8]])

要为ndarray情况沿每个轴指定不同的仓大小,可以沿每个轴迭代使用np.add.reduceat,如so-

def bucket(x, bin_size):
    ndims = x.ndim
    out = x.copy()
    for i in range(ndims):
        idx = np.append(0,np.cumsum(bin_size[i][:-1]))
        out = np.add.reduceat(out,idx,axis=i)
    return out

样品运行-

In [126]: x
Out[126]: 
array([[165, 107, 133,  82, 199],
       [ 35, 138,  91, 100, 207],
       [ 75,  99,  40, 240, 208],
       [166, 171,  78,   7, 141]])
In [127]: bucket(x, bin_size = [[2, 2],[3, 2]])
Out[127]: 
array([[669, 588],
       [629, 596]])
#  [2, 2] are the bin sizes along axis=0
#  [3, 2] are the bin sizes along axis=1
# array([[165, 107, 133, | 82, 199],
#        [ 35, 138,  91, | 100, 207],
# -------------------------------------
#        [ 75,  99, 40,  | 240, 208],
#        [166, 171, 78,  | 7, 141]])
In [128]: x[:2,:3].sum()
Out[128]: 669
In [129]: x[:2,3:].sum()
Out[129]: 588
In [130]: x[2:,:3].sum()
Out[130]: 629
In [131]: x[2:,3:].sum()
Out[131]: 596

最新更新