使用 NumPy 更优化的数学函数

Question

对于我的班级，我需要使用 NumPy 编写更优化的数学函数。问题是，当使用 NumPy 时，我的解决方案在使用原生 Python 时会更慢。

1. 对数组的所有元素求和并求和的函数

Python：

def cube(x):
    result = 0
    for i in range(len(x)):
        result += x[i] ** 3
    return result

我的，使用 NumPy（慢 15-30%）：

def cube(x):
    it = numpy.nditer([x, None])
    for a, b in it:
        b[...] = a*a*a
    return numpy.sum(it.operands[1])

1. 一些随机计算函数

Python：

def calc(x):
    m = sum(x) / len(x)
    result = 0

    for i in range(len(x)):
        result += (x[i] - m)**4

    return result / len(x)

NumPy（>10 倍慢）：

def calc(x):
    m = numpy.mean(x)
    result = 0
    for i in range(len(x)):
        result += numpy.power((x[i] - m), 4)
    return result / len(x)

我不知道如何解决这个问题，到目前为止我已经尝试过 NumPy 的随机函数

Answer 1

详细说明 cmets 中所说的内容：

Numpy 的强大之处在于能够在快速 c/fortran 中完成所有循环，而不是慢速 Python 循环。例如，如果您有一个数组x，并且您想计算该数组中每个值的平方，您可以这样做

y = []
for value in x:
    y.append(value**2)

甚至（使用列表理解）

y = [value**2 for value in x]

但是如果你可以在 numpy 中进行所有循环，它会更快

y = x**2

（假设 x 已经是一个 numpy 数组）。

因此，对于您的示例，在 numpy 中执行此操作的正确方法是

1.

def sum_of_cubes(x):
    result = 0
    for i in range(len(x)):
        result += x[i] ** 3
    return result

def sum_of_cubes_numpy(x):
    return (x**3).sum()

def calc(x):
    m = sum(x) / len(x)
    result = 0

    for i in range(len(x)):
        result += (x[i] - m)**4

    return result / len(x)

def calc_numpy(x):
    m = numpy.mean(x)  # or just x.mean()
    return numpy.sum((x - m)**4) / len(x)

请注意，我假设输入 x 已经是一个 numpy 数组，而不是常规 Python 列表：如果您有一个列表 lst，您可以使用 arr = numpy.array(lst) 从它创建一个数组。

Answer 2

In [337]: def cube(x):
     ...:     result = 0
     ...:     for i in range(len(x)):
     ...:         result += x[i] ** 3
     ...:     return result
     ...: 

nditer 不是一个好的 numpy 迭代器，至少在 Python 级代码中使用时不是这样。它实际上只是编写编译代码的垫脚石。它的文档需要更好的免责声明。

In [338]: def cube1(x):
     ...:     it = numpy.nditer([x, None])
     ...:     for a, b in it:
     ...:         b[...] = a*a*a
     ...:     return numpy.sum(it.operands[1])
     ...: 
In [339]: cube(list(range(10)))
Out[339]: 2025
In [340]: cube1(list(range(10)))
Out[340]: 2025
In [341]: cube1(np.arange(10))
Out[341]: 2025

更直接的numpy 迭代：

In [342]: def cube2(x):
     ...:     it = [a*a*a for a in x]
     ...:     return numpy.sum(it)
     ...: 

更好的全数组代码。正如 sum 可以对整个数组起作用一样，power 也适用于整个数组。

In [343]: def cube3(x):
     ...:     return numpy.sum(x**3)
     ...: 
In [344]: cube2(np.arange(10))
Out[344]: 2025
In [345]: cube3(np.arange(10))
Out[345]: 2025

做一些计时：

列表参考：

In [346]: timeit cube(list(range(1000)))
438 µs ± 9.87 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

慢的nditer：

In [348]: timeit cube1(np.arange(1000))
2.8 ms ± 5.65 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

部分numpy：

In [349]: timeit cube2(np.arange(1000))
520 µs ± 20 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

我可以通过传递一个列表而不是一个数组来改进它的时间。列表的迭代速度更快。

In [352]: timeit cube2(list(range(1000)))
229 µs ± 9.53 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

但是“纯”numpy 版本的时间让所有这些都失去了意义：

In [350]: timeit cube3(np.arange(1000))
23.6 µs ± 128 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

一般规则是应用于numpy 数组的numpy 方法最快。但如果必须循环，通常最好使用列表。

有时纯numpy 方法会创建非常大的临时数组。那么内存管理的复杂性会降低性能。在这种情况下，对复杂任务进行适度的迭代可能是最好的。