Numpy：初学者

Question

我正在尝试学习 nditer 以可能用于加快我的应用程序。在这里，我尝试制作一个有趣的重塑程序，该程序将采用大小为 20 的数组并将其重塑为 5x4 数组：

myArray = np.arange(20)
def fi_by_fo_100(array):
    offset = np.array([0, 4, 8, 12, 16])
    it = np.nditer([offset, None],
                      flags=['reduce_ok'],
                      op_flags=[['readonly'],
                                ['readwrite','allocate']],
                      op_axes=[None, [0,1,-1]],
                      itershape=(-1, 4, offset.size))

    while not it.finished:
        indices = np.arange(it[0],(it[0]+4), dtype=int)
        info = array.take(indices)
        '''Just for fun, we'll perform an operation on data.\
           Let's shift it to 100'''
        info = info + 81
        it.operands[1][...]=info
        it.iternext()
    return it.operands[1]

test = fi_by_fo_100(myArray)
>>> test
array([[ 97,  98,  99, 100]])

显然，程序将每个结果都覆盖到一行中。所以我尝试使用nditer的索引功能，但仍然没有骰子。

flags=['reduce_ok','c_iter'] --> it.operands[1][it.index][...]=info =
IndexError: index out of bounds

flags=['reduce_ok','c_iter'] --> it.operands[1][it.iterindex][...]=info =
IndexError: index out of bounds

flags=['reduce_ok','multi_iter'] --> it.operands[1][it.multi_index][...]=info =
IndexError: index out of bounds

it[0][it.multi_index[1]][...]=info =
IndexError: 0-d arrays can't be indexed

...等等。我错过了什么？提前致谢。

奖金问题

我刚刚遇到this nice article on nditer。我可能是 Numpy 的新手，但这是我第一次看到 Numpy 速度基准如此落后。我的理解是人们选择 Numpy 是因为它的数值速度和能力，但迭代是其中的一部分，不是吗？如果它这么慢，nditer 有什么意义？

Answer 1

打印出一路上发生的事情确实有助于分解事情。

首先，让我们用这个替换整个循环：

i = 0
while not it.finished:
    i += 1
print i

它将打印 20，而不是 5。那是因为您正在进行 5x4 迭代，而不是 5x1。

那么，为什么这甚至接近工作？好吧，让我们更仔细地看一下循环：

while not it.finished:
    print '>', it.operands[0], it[0]
    indices = np.arange(it[0],(it[0]+4), dtype=int)
    info = array.take(indices)
    info = info + 81
    it.operands[1][...]=info
    print '<', it.operands[1], it[1]

您会看到前五个循环经过[0 4 8 12 16] 五次，生成[[81 82 83 84]]，然后是[[85 86 87 88]]，等等。然后接下来的五个循环做同样的事情，一次又一次。

这也是您的 c_index 解决方案不起作用的原因 - 因为 it.index 的范围是 0 到 19，而您在 it.operands[1] 中没有 20 个。

如果您正确执行了 multi_index 并忽略了列，则可以完成这项工作……但是，您仍然会进行 5x4 迭代，只是将每个步骤重复 4 次，而不是进行您想要的 5x1 迭代。

您的it.operands[1][...]=info 每次循环都会用 5x1 行替换整个输出。一般来说，你不应该对it.operands[1] 做任何事情——nditer 的全部意义在于你只需处理好每个it[1]，最后的it.operands[1] 就是结果。

当然，对行进行 5x4 迭代是没有意义的。对单个值进行 5x4 迭代，或者对行进行 5x1 迭代。

如果你想要前者，最简单的方法是重塑输入数组，然后迭代：

it = np.nditer([array.reshape(5, -1), None],
               op_flags=[['readonly'],
                         ['readwrite','allocate']])
for a, b in it:
    b[...] = a + 81
return it.operands[1]

但这当然很愚蠢——它只是一种更慢、更复杂的写作方式：

return array+81

建议“编写自己的reshape 的方法是先调用reshape，然后……”会有点愚蠢。

那么，您想遍历行，对吗？

让我们通过去掉allocate 并明确创建一个 5x4 数组来简化一下：

outarray = np.zeros((5,4), dtype=array.dtype)
offset = np.array([0, 4, 8, 12, 16])
it = np.nditer([offset, outarray],
               flags=['reduce_ok'],
               op_flags=[['readonly'],
                         ['readwrite']],
               op_axes=[None, [0]],
               itershape=[5])

while not it.finished:
    indices = np.arange(it[0],(it[0]+4), dtype=int)
    info = array.take(indices)
    '''Just for fun, we'll perform an operation on data.\
       Let's shift it to 100'''
    info = info + 81
    it.operands[1][it.index][...]=info
    it.iternext()
return it.operands[1]

这有点滥用nditer，但至少它做了正确的事。

由于您只是对源代码进行 1D 迭代而基本上忽略了第二次迭代，因此在此处使用 nditer 确实没有充分的理由。如果您需要对多个数组进行锁步迭代，for a, b in nditer([x, y], …) 比迭代x 并使用索引访问y 更干净——就像for a, b in zip(x, y) 在numpy 之外一样。如果您需要遍历多维数组，nditer 通常比其他方法更干净。但是在这里，您真正要做的只是遍历[0, 4, 8, 16, 20]，对结果进行处理，然后将其复制到另一个array。

另外，正如我在 cmets 中提到的，如果您发现自己在 numpy 中使用迭代，那么您通常做错了什么。 numpy 的所有速度优势都来自于让它在本机 C/Fortran 或更低级别的向量操作中执行紧密循环。一旦你在 arrays 上循环，你实际上只是在用稍微好一点的语法做慢速 Python 数字：

import numpy as np
import timeit

def add10_numpy(array):
    return array + 10

def add10_nditer(array):
    it = np.nditer([array, None], [],
                   [['readonly'], ['writeonly', 'allocate']])
    for a, b in it:
        np.add(a, 10, b)
    return it.operands[1]

def add10_py(array):
    x, y = array.shape
    outarray = array.copy()
    for i in xrange(x):
        for j in xrange(y):
            outarray[i, j] = array[i, j] + 10
    return out array

myArray = np.arange(100000).reshape(250,-1)

for f in add10_numpy, add10_nditer, add10_py:
    print '%12s: %s' % (f.__name__, timeit.timeit(lambda: f(myArray), number=1))

在我的系统上，打印：

 add10_numpy: 0.000458002090454
add10_nditer: 0.292730093002
    add10_py: 0.127345085144

这显示了不必要地使用nditer 的成本。