1

我在 Python 中有这段代码

for i in range(len(ax)):
  for j in range(len(rx)):
    x = ax[i] + rx[j]
    y = ay[i] + ry[j]
    A[x,y] = A[x,y] + 1

在哪里

A.shape = (N,M)
ax.shape = ay.shape = (L)
rx.shape = ry.shape = (K)

我想矢量化或以其他方式使其更高效，即更快，并且如果可能的话在内存消耗上更经济。这里，我的 ax 和 ay 指的是数组 A 的绝对元素，而 rx 和 ay 是相对坐标。所以，我正在更新计数器数组 A。

我的表 A 可以是 1000x1000，而 ax,ay 是 100x1，cx,cy 是 300x1。整个事情都在循环中，最好是优化的代码不会继续创建 A 大小的大表。

这个问题与我之前提出的问题有关，但由于增量的工作方式，它并不直接适用于这种情况。这是一个例子。

这段代码正是我想要的：

import numpy as np
A = np.zeros((4,5))
ax = np.arange(1,3)
ay = np.array([1,1])
rx = np.array([-1,0,0])
ry = np.array([0,0,0])

for i in range(len(ax)):
    for j  in range(len(rx)):
        x = ax[i] + rx[j]
        y = ay[i] + ry[j]
        print(x,y)
        A[x,y] = A[x,y] + 1
A
array([[ 0.,  1.,  0.,  0.,  0.],
       [ 0.,  3.,  0.,  0.,  0.],
       [ 0.,  2.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.]])

但是，下面的代码不起作用，因为当我们递增一个数组时，它会预先计算数组的右侧：

import numpy as np
A = np.zeros((4,5))
ax = np.arange(1,3)
ay = np.array([1,1])
rx = np.array([-1,0])
ry = np.array([0,0])

x = ax + rx[:,np.newaxis]
y = ay + ry[:,np.newaxis]
A[x,y] = A[x,y] + 1

A

array([[ 0.,  1.,  0.,  0.,  0.],
       [ 0.,  1.,  0.,  0.,  0.],
       [ 0.,  1.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.]])

这个解决方案在数字的正确性方面起作用，但它不是最快的，可能是因为 np.add.at() 函数没有缓冲：

import numpy as np
A = np.zeros((4,5))
ax = np.arange(1,3)
ay = np.array([1,1])
rx = np.array([-1,0,0])
ry = np.array([0,0,0])

x = ax + rx[:,np.newaxis]
y = ay + ry[:,np.newaxis]
np.add.at(A,[x,y],1)

A
4

1 回答 1

0

这是一个利用broadcasting，获得线性索引，然后将其馈送到非常有效np.bincount的合并求和 -

m,n = 4,5 # shape of output array
X = ax[:,None] + rx
Y = ay[:,None] + ry
Aout = np.bincount((X*n + Y).ravel(), minlength=m*n).reshape(m,n)

替代方案之一np.flatnonzero-

idx = (X*n + Y).ravel()
idx.sort()
m = np.r_[True,idx[1:] != idx[:-1],True]
A.ravel()[idx[m[:-1]]] = np.diff(np.flatnonzero(m))

如果您要A迭代添加，请在最后一步替换=为+=那里。

于 2017-12-22T11:07:10.767 回答