这是正常和预期的行为。应用“避免使用 numpy 的循环”语句太简单了。如果您正在处理内部循环,那么它(几乎)总是正确的。但是在外部循环的情况下(就像你的情况一样),还有更多的例外。特别是如果替代方案是使用广播,因为这会通过使用更多内存来加速您的操作。
只是为“避免使用 numpy 循环”语句添加一些背景知识:
NumPy 数组存储为具有c类型的连续数组。Pythonint
与 C 不同int
!因此,每当您遍历数组中的每个项目时,您都需要从数组中插入该项目,将其转换为 Python int
,然后对它进行任何您想做的事情,最后您可能需要再次将其转换为 ac 整数(称为装箱并拆箱值)。例如,您想sum
使用 Python 对数组中的项目:
import numpy as np
arr = np.arange(1000)
%%timeit
acc = 0
for item in arr:
acc += item
# 1000 loops, best of 3: 478 µs per loop
你最好使用numpy:
%timeit np.sum(arr)
# 10000 loops, best of 3: 24.2 µs per loop
即使您将循环推送到 Python C 代码中,您也远离 numpy 性能:
%timeit sum(arr)
# 1000 loops, best of 3: 387 µs per loop
此规则可能有例外,但这些将非常稀少。至少只要有一些等效的 numpy 功能。因此,如果您要遍历单个元素,那么您应该使用 numpy.
有时一个普通的 python 循环就足够了。它没有被广泛宣传,但与 Python 函数相比,numpy 函数具有巨大的开销。例如考虑一个 3 元素数组:
arr = np.arange(3)
%timeit np.sum(arr)
%timeit sum(arr)
哪个会更快?
解决方案:Python函数比numpy解决方案执行得更好:
# 10000 loops, best of 3: 21.9 µs per loop <- numpy
# 100000 loops, best of 3: 6.27 µs per loop <- python
但这与您的示例有什么关系?事实上并不是那么多,因为你总是在数组上使用 numpy 函数(不是单个元素,甚至不是几个元素),所以你的内部循环已经使用了优化的函数。这就是为什么两者的性能大致相同(+/- 10 倍,元素很少,到 2 倍,大约 500 个元素)。但这并不是真正的循环开销,而是函数调用开销!
您的循环解决方案
使用line-profiler和 a resolution = 100
:
def fun_func(tim, prec, values):
for i, ti in enumerate(tim):
values[i] = np.sum(np.sin(prec * ti))
%lprun -f fun_func fun_func(tim, prec, values)
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 101 752 7.4 5.7 for i, ti in enumerate(tim):
3 100 12449 124.5 94.3 values[i] = np.sum(np.sin(prec * ti))
95% 都花在循环内,我什至将循环体分成几个部分来验证这一点:
def fun_func(tim, prec, values):
for i, ti in enumerate(tim):
x = prec * ti
x = np.sin(x)
x = np.sum(x)
values[i] = x
%lprun -f fun_func fun_func(tim, prec, values)
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 101 609 6.0 3.5 for i, ti in enumerate(tim):
3 100 4521 45.2 26.3 x = prec * ti
4 100 4646 46.5 27.0 x = np.sin(x)
5 100 6731 67.3 39.1 x = np.sum(x)
6 100 714 7.1 4.1 values[i] = x
消费者在np.multiply
, np.sin
,的时间np.sum
在这里,您可以通过比较他们每次调用的时间与他们的开销来轻松检查:
arr = np.ones(1, float)
%timeit np.sum(arr)
# 10000 loops, best of 3: 22.6 µs per loop
因此,与计算运行时相比,只要计算函数调用开销很小,您就会有类似的运行时。即使有 100 个项目,您也非常接近开销时间。诀窍是知道他们在什么时候收支平衡。对于 1000 个项目,调用开销仍然很大:
%lprun -f fun_func fun_func(tim, prec, values)
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 1001 5864 5.9 2.4 for i, ti in enumerate(tim):
3 1000 42817 42.8 17.2 x = prec * ti
4 1000 119327 119.3 48.0 x = np.sin(x)
5 1000 73313 73.3 29.5 x = np.sum(x)
6 1000 7287 7.3 2.9 values[i] = x
但是与resolution = 5000
运行时相比,开销非常低:
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 5001 29412 5.9 0.9 for i, ti in enumerate(tim):
3 5000 388827 77.8 11.6 x = prec * ti
4 5000 2442460 488.5 73.2 x = np.sin(x)
5 5000 441337 88.3 13.2 x = np.sum(x)
6 5000 36187 7.2 1.1 values[i] = x
当您在每次通话中花费 500 美元np.sin
时,您不再关心 20 美元的开销。
需要注意的是:line_profiler
可能包括每行的一些额外开销,也可能包括每个函数调用,因此函数调用开销变得可以忽略不计的点可能会更低!!!
您的广播解决方案
我从分析第一个解决方案开始,让我们对第二个解决方案做同样的事情:
def fun_func(tim, prec, values):
x = tim[:, np.newaxis]
x = x * prec
x = np.sin(x)
x = np.sum(x, axis=1)
return x
再次使用 line_profiler resolution=100
:
%lprun -f fun_func fun_func(tim, prec, values)
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 1 27 27.0 0.5 x = tim[:, np.newaxis]
3 1 638 638.0 12.9 x = x * prec
4 1 3963 3963.0 79.9 x = np.sin(x)
5 1 326 326.0 6.6 x = np.sum(x, axis=1)
6 1 4 4.0 0.1 return x
这已经大大超过了开销时间,因此我们最终比循环快了 10 倍。
我还对以下内容进行了分析resolution=1000
:
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 1 28 28.0 0.0 x = tim[:, np.newaxis]
3 1 17716 17716.0 14.6 x = x * prec
4 1 91174 91174.0 75.3 x = np.sin(x)
5 1 12140 12140.0 10.0 x = np.sum(x, axis=1)
6 1 10 10.0 0.0 return x
并与precision=5000
:
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def fun_func(tim, prec, values):
2 1 34 34.0 0.0 x = tim[:, np.newaxis]
3 1 333685 333685.0 11.1 x = x * prec
4 1 2391812 2391812.0 79.6 x = np.sin(x)
5 1 280832 280832.0 9.3 x = np.sum(x, axis=1)
6 1 14 14.0 0.0 return x
1000 大小仍然更快,但正如我们在那里看到的那样,在循环解决方案中调用开销仍然不可忽略。但是resolution = 5000
每个步骤所花费的时间几乎相同(有些慢一些,有些快但总体上非常相似)
另一个影响是,当你做乘法时,实际的广播变得很重要。即使使用非常智能的 numpy 解决方案,这仍然包括一些额外的计算。如resolution=10000
您所见,广播乘法开始占用与循环解决方案相关的更多“% 时间”:
Line # Hits Time Per Hit % Time Line Contents
==============================================================
1 def broadcast_solution(tim, prec, values):
2 1 37 37.0 0.0 x = tim[:, np.newaxis]
3 1 1783345 1783345.0 13.9 x = x * prec
4 1 9879333 9879333.0 77.1 x = np.sin(x)
5 1 1153789 1153789.0 9.0 x = np.sum(x, axis=1)
6 1 11 11.0 0.0 return x
Line # Hits Time Per Hit % Time Line Contents
==============================================================
8 def loop_solution(tim, prec, values):
9 10001 62502 6.2 0.5 for i, ti in enumerate(tim):
10 10000 1287698 128.8 10.5 x = prec * ti
11 10000 9758633 975.9 79.7 x = np.sin(x)
12 10000 1058995 105.9 8.6 x = np.sum(x)
13 10000 75760 7.6 0.6 values[i] = x
但是除了实际花费的时间之外还有另一件事:内存消耗。您的循环解决方案需要O(n)
内存,因为您总是处理n
元素。然而,广播解决方案需要O(n*n)
内存。如果你使用你的循环,你可能需要等待一段时间,resolution=20000
但它仍然只需要8bytes/element * 20000 element ~= 160kB
,但你需要广播~3GB
。这忽略了常数因素(如临时数组又名中间数组)!假设你走得更远,你会很快耗尽内存!
是时候再次总结一下要点了:
- 如果您对 numpy 数组中的单个项目执行 python 循环,那么您做错了。
- 如果您遍历 numpy-array 的子数组,请确保每个循环中的函数调用开销与函数中花费的时间相比可以忽略不计。
- 如果您广播 numpy 数组,请确保您不会耗尽内存。
但是关于优化最重要的一点仍然是:
最后一个想法:
如果numpy或scipy中没有现有的解决方案,则可以使用cython、numba或numexpr轻松实现此类需要循环或广播的功能。
例如,将循环解决方案的内存效率与低广播解决方案的速度相结合的 numba 函数resolutions
如下所示:
from numba import njit
import math
@njit
def numba_solution(tim, prec, values):
size = tim.size
for i in range(size):
ti = tim[i]
x = 0
for j in range(size):
x += math.sin(prec[j] * ti)
values[i] = x
正如评论中所指出的,还numexpr
可以非常快速地评估广播计算并且不需要O(n*n)
内存:
>>> import numexpr
>>> tim_2d = tim[:, np.newaxis]
>>> numexpr.evaluate('sum(sin(tim_2d * prec), axis=1)')