我可以通过numpy避免Python循环开销的动态编程吗？

Question

我需要以下问题的Pythonic循环开销的帮助：我正在写一个函数来计算像素流算法，它是2D Numpy数组上的经典动态编程算法。它要求：

1）访问数组中的所有元素至少一次这样的：

for x in range(xsize): 
    for y in range(ysize): 
     updateDistance(x,y) 

2）潜在的下列基于它看起来像一个元素的邻居的值的元素的一个路径这

while len(workingList) > 0: 
    x,y = workingList.pop() 
    #if any neighbors of x,y need calculation, push x,y and neighbors on workingList 
    #else, calculate flow on pixels as a sum of flow on neighboring pixels 

不幸的是，我似乎得到了很多Python的循环开销对＃1，即使调用updateDistance是pass。我认为这是一个足够经典的算法，在Python中必须有一个很好的方法来避免一些循环开销。我也担心，如果我能修好＃1，我会在＃2的循环中被淘汰。

任何关于快速遍历2D numpy数组中的元素并可能更新该数组中的元素链的建议？

编辑：冲洗出来的＃2

更多的细节看来我也许能向量化​​的第一环，或许与矢量化的np.meshgrid电话。

该循环部分有点复杂，但这里是一个简化版本。我很关心循环和索引到相邻元素：

#A is a 2d cost matrix 
workingList = [(x,y)] 
while len(workingList) > 0: 
    x,y = workingList.pop() 
    neighborsToCalculate = [] 
    for n in neighborsThatNeedCalculation(x,y): #indexes A to check neighbors of (x,y) 
     neighborsToCalculate.append(n) 
    if len(neighborstToCalculate) != 0: 
     workingList.append((x,y)) 
     workingList.extend(neighborsToCalculate) 
    else: 
     for xn,yn in neighbors(x,y): 
      A[x,y] += 1+A[xn,yn] 

这是一个经典的宽度优先搜索问题。如果这可以并行，那将是非常好的。它可能不能以其目前的形式，因为它遵循了一条道路，但我很乐意提供建议。

Answer 1

如果你在算法中使用python循环，你将不会从numpy获得任何速度提升。你需要并行处理你的问题。

在图像处理中，并行化意味着在所有像素上使用相同的函数，即使用内核。在numpy的，而不是这样做的：

for x in range(xsize): 
    for y in range(ysize): 
     img1[y, x] = img2[y, x] + img3[y, x] 

你这样做：

img1 = img2 + img3 # add 2 images pixelwise 

使回路发生在C。事实上，你有一个每个像素长度未知的邻居列表，这使得你的问题很难以这种方式进行并行化。你应该修改你的问题（你可以更具体一点关于你的算法吗？），或者使用其他语言，比如cython。

编辑：

你不会得到好处numpy的不改变你的算法。 Numpy允许你执行线性代数运算。执行任意操作时，您无法避免使用此库的循环开销。

优化这个，你可以考虑：

• 切换到喜欢用Cython另一种语言（这是一家专业从事Python扩展）摆脱循环成本

• 优化算法：如果您可以使用线性代数运算（这取决于函数neighborsThatNeedCalculation）获得相同的结果，但可以使用numpy，但需要制定新的体系结构。

• 使用像MapReduce这样的并行化技术。有了python，你可以使用一组工作者（可用于多处理模块），如果你切换到另一种语言，你将获得更多的速度收益，因为python将有其他瓶颈。

在这种情况下你想要的东西易于安装和集成，而你只需要拥有C-般的表演，我强烈建议用Cython如果你不能返工你的算法。

Answer 2

对于第一部分，您可以使用numpy.vectorize，但只有在无法使用数组操作来实现updateDistance的功能时才会这样做。这里有一个例子：

import numpy as np  
updateDistance = np.vectorize(lambda x: x + 1) # my updateDistance increments 

在现实中，如果这是你正在尝试做的工作，只是做a + 1。因此，如果我们采取的一些的阵列和应用updateDistance：

>>> a = np.ones((3,3)) 
>>> updateDistance(a) 
array([[ 2., 2., 2.], 
     [ 2., 2., 2.], 
     [ 2., 2., 2.]]) 

至于第二部分，我不认为我了解细节不够好提出一个更好的选择。这听起来像是你需要反复观察最近的邻居，所以我怀疑你至少可以改进if-else中的东西。

---

更新：时机与第一部分。

注意：这些时间是在我的机器上完成的，没有试图规范化环境。

python -mtimeit 'import numpy as np' 'n = 100' 'a = np.ones((n, n))' 'b = np.zeros((n, n))' 'for x in range(n): ' ' for y in range(n):' '  b[x,y] = a[x,y] + 1' 

的np.vectorize时间与生成：

循环时间与生成

python -mtimeit 'import numpy as np' 'n = 100' 'a = np.ones((n, n))' 'updateDistance = np.vectorize(lambda x: x + 1)' 'b = updateDistance(a)' 

在两种情况下，n = 100导致100×100阵列。根据需要替换100。

Array size Loop version np.vectorize version np.vectorize speed up 
100 x 100  20.2 msec  2.6 msec    7.77x 
200 x 200  81.8 msec  10.4 msec    7.87x 
400 x 400  325 msec  42.6 msec    7.63x 

最后，在np.vectorize例如与单纯使用数组操作，你可以做比较：

python -mtimeit 'import numpy as np' 'n = 100' 'a = np.ones((n, n))' 'a += 1' 

在我的机器，这产生了以下内容。

Array size Array operation version Speed up over np.vectorize version 
100 x 100  23.6 usec     110.2x 
200 x 200  79.7 usec     130.5x 
400 x 400  286 usec     149.0x 

总之，在使用np.vectorize，而不是循环的优点，但有一个更大的使用数组操作，如果可以实现的updateDistance功能诱因。

Answer 3

你应该考虑使用C-extension/Cython。 如果你留在Python中，一个主要的改进可以通过更换来实现：

for xn,yn in neighbors(x,y): 
     A[x,y] += 1+A[xn,yn] 

有：

n = neighbors(x,y) 
A[x,y] += len(n)+sum(A[n]) 

neighbors应该返回指标，不标。