矢量化调用数组第三维的numpy函数



我有一个3D数字数组data,其中维度ab表示图像的分辨率,c是图像/帧编号。我想在c维度上的每个像素(ab组合(上调用np.histogram,并使用维度为(a, b, BINS)的输出数组。我已经用嵌套循环完成了这项任务,但如何将此操作向量化呢?

hists = np.zeros((a, b, BINS))
for row in range(a):
for column in range(b):
hists[row, column, :] = np.histogram(data[row, column, :], bins=BINS)[0]

我相信这个解决方案是微不足道的,尽管如此,我们还是感谢所有的帮助:(

np.histogram在扁平阵列上进行计算。但是,您可以使用np.apply_along_axis

np.apply_along_axis(lambda a: np.histogram(a, bins=BINS)[0], 2, data)

这是一个有趣的问题。

制作最小工作示例(MWE(

这应该是在SO上提问的主要习惯。

a, b, c = 2, 3, 4
data = np.random.randint(10, size=(a, b, c))
hists = np.zeros((a, b, c), dtype=int)
for row in range(a):
for column in range(b):
hists[row, column, :] = np.histogram(data[row, column, :], bins=c)[0]
data
>>> array([[[6, 4, 3, 3],
[7, 3, 8, 0],
[1, 5, 8, 0]],
[[5, 5, 7, 8],
[3, 2, 7, 8],
[6, 8, 8, 0]]])
hists
>>> array([[[2, 1, 0, 1],
[1, 1, 0, 2],
[2, 0, 1, 1]],
[[2, 0, 1, 1],
[2, 0, 0, 2],
[1, 0, 0, 3]]])

尽可能简单(但仍然有效(

你可以消除一个循环并简化它:

new_data = data.reshape(a*b, c)
new_hists = np.zeros((a*b, c), dtype=int)
for row in range(a*b):
new_hists[row, :] = np.histogram(new_data[row, :], bins=c)[0]
new_hists
>>> array([[2, 1, 0, 1],
[1, 1, 0, 2],
[2, 0, 1, 1],
[2, 0, 1, 1],
[2, 0, 0, 2],
[1, 0, 0, 3]])
new_data
>>> array([[6, 4, 3, 3],
[7, 3, 8, 0],
[1, 5, 8, 0],
[5, 5, 7, 8],
[3, 2, 7, 8],
[6, 8, 8, 0]])

你能找到类似的问题并使用它们解决方案的关键点吗

一般来说,你不能对循环中正在进行的事情进行矢量化:

for row in array:
some_operation(row)

除了这种情况,您可以对扁平阵列调用另一个矢量化操作,然后将其移回初始形状:

arr = array.ravel()
another_operation(arr)
out = arr.reshape(array.shape)

看起来你很幸运有了np.histogram,因为我敢肯定以前也做过类似的事情。

最终解决方案

new_data = data.reshape(a*b, c)
m, M = new_data.min(axis=1), new_data.max(axis=1)
bins = (c * (new_data - m[:, None]) // (M-m)[:, None])
out = np.zeros((a*b, c+1), dtype=int)
advanced_indexing = np.repeat(np.arange(a*b), c), bins.ravel()
np.add.at(out, advanced_indexing, 1)
out.reshape((a, b, -1))
>>> array([[[2, 1, 0, 0, 1],
[1, 1, 0, 1, 1],
[2, 0, 1, 0, 1]],
[[2, 0, 1, 0, 1],
[2, 0, 0, 1, 1],
[1, 0, 0, 1, 2]]])

请注意,它在每个直方图中添加了一个额外的bin,并在其中放置最大值,但如果需要,我希望它不难修复。

最新更新