如何有效地操作一个大的numpy数组

Posted 2023-03-12

【中文标题】如何有效地操作一个大的numpy数组

【英文标题】：How to efficiently operate a large numpy array

【发布时间】：2019-03-31 15:53:39

【问题描述】：

我有一段代码，它基于一个大的 numpy 数组，然后操作另一个数组。因为这是一个非常大的数组，请问是否有一种有效的方法可以实现我的目标？ （我认为高效的方式应该是直接对数组进行操作，而不是通过for循环）。

提前致谢，请在下面找到我的代码：

N = 1000000000
rand = np.random.rand(N)
beta = np.zeros(N)
for i in range(0, N):
    if rand[i] < 0.5:
        beta[i] = 2.0*rand[i]
    else:
        beta[i] = 1.0/(2.0*(1.0-rand[i]))

【问题讨论】：

你可以试试scipy的稀疏矩阵。查看此post。

对于这么大的矩阵，这可能仍然太慢，但可能更好：np.where(rand < 0.5, 2.0*rand, 1.0/(2.0*(1.0-rand)))

@EricWang，我看不出稀疏矩阵与这个问题的相关性rand 没有高比例的零，beta 也没有。

一种更快的方法，几乎​​不需要更改您的代码（我不会初始化一个完全用零覆盖并添加 @nb.njit(error_model="numpy) 的数组）将是 ***.com/q/52046301/4045774 .

【参考方案1】：

通过在 Python 中执行处理，您在这里基本上失去了 numpy 的效率。 numpy 的想法是处理 bulk 中的项目，因为它在 C++ 中具有执行实际处理的效率算法。您可以将 numpy 的 Python 结尾更多地视为“接口”。

现在回答你的问题，我们基本上可以先构造一个介于 0 和 2 之间的随机数数组，方法是将它与 2 相乘：

rand = 2.0 * np.random.rand(N)

接下来我们可以使用np.where(..) [numpy-doc]，它就像一个条件选择器：我们在这里传递三个“数组”：第一个是一个布尔数组，它编码“条件”的真实性，第二个是数组在相关条件为真时要填写的值，第三个值是在条件为假时要插入的值的数组，所以我们可以这样写：

N = 1000000000
rand = 2 * np.random.rand(N)
beta = np.where(rand < 1.0, rand, 1.0 / (2.0 - rand))

【讨论】：

非常感谢，根据答案，有两种方法：1）在我定义一个函数后使用 numpy 的矢量化，2）使用 where。我在我的电脑上用 N=1x10^8 进行了测试，第一种方法需要 30 秒，第二种方法需要 1.3 秒。如果我没有先乘以 2，则需要 1.9 秒。非常感谢。

【参考方案2】：

N = 1000000000 为我造成了MemoryError。以最小的例子减少到 100。 您可以使用np.where routine。

在这两种情况下，从根本上讲，您都是在遍历数组并应用一个函数。 然而，np.where 使用了一种更快的循环（它基本上是编译代码），而你的“python”循环被解释，因此对于一个大 N 来说真的很慢。

这是一个实现示例。

N = 100
rand = np.random.rand(N)
beta = np.where(rand < 0.5,  2.0 * rand, 1.0/(2.0*(1.0-rand))

【参考方案3】：

正如其他答案所指出的，在 Python 循环中迭代 numpy 数组的元素几乎总是应该（并且可以）避免。在大多数情况下，从 Python 循环到数组操作的速度提升约为 100 倍。

但是，如果性能绝对关键，您通常可以使用 Cython 挤出 2 倍到 10 倍之间的另一个因素（根据我的经验）。 这是一个例子：

%%cython
cimport numpy as np
import numpy as np
cimport cython
from cython cimport floating

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.cdivision(True)
cpdef np.ndarray[floating, ndim=1] beta(np.ndarray[floating, ndim=1] arr):
    cdef:
        Py_ssize_t i
        Py_ssize_t N = arr.shape[0]
        np.ndarray[floating, ndim=1] result = np.zeros(N)

    for i in range(N):
        if arr[i] < 0.5:
            result[i] = 2.0*arr[i]
        else:
            result[i] = 1.0/(2.0*(1.0-arr[i]))

    return result

然后您可以将其称为beta(rand)。 如您所见，这允许您使用原始循环结构，但现在使用高效的类型化本机代码。与 np.where 相比，我得到了约 2.5 倍的加速。

需要注意的是，在许多情况下，与 numpy 中的单行代码相比，这不值得付出额外的努力——但这很可能是性能至关重要的地方。

【讨论】：

非常感谢，这也是我不知道Numpy高效使用时的一个好答案。