将函数映射到 numpy 数组的最有效方法是什么?我目前正在做:
import numpy as np
x = np.array([1, 2, 3, 4, 5])
# Obtain array of square of each element in x
squarer = lambda t: t ** 2
squares = np.array([squarer(xi) for xi in x])
但是,这可能效率很低,因为我使用列表推导将新数组构造为 Python 列表,然后再将其转换回 numpy 数组。我们能做得更好吗?
squarer(x)
怎么样?
squarer(x)
将对数组元素应用 squarer
函数并返回一个数组,其中包含单个 squarer(element)
调用的结果。我写这个是因为“只有 squarer(x) 怎么样?”乍一看还不够清楚。
我已经用 perfplot
(我的一个小项目)测试了所有建议的方法以及 np.array(list(map(f, x)))
。
消息 #1:如果您可以使用 numpy 的本机函数,请执行此操作。
如果您尝试矢量化的函数已经 被矢量化(如原始帖子中的 x**2
示例),则使用它比其他任何方法都快(请注意对数刻度):
https://i.stack.imgur.com/IDYsC.png
如果您确实需要矢量化,那么您使用哪种变体并不重要。
https://i.stack.imgur.com/V30x6.png
重现绘图的代码:
import numpy as np
import perfplot
import math
def f(x):
# return math.sqrt(x)
return np.sqrt(x)
vf = np.vectorize(f)
def array_for(x):
return np.array([f(xi) for xi in x])
def array_map(x):
return np.array(list(map(f, x)))
def fromiter(x):
return np.fromiter((f(xi) for xi in x), x.dtype)
def vectorize(x):
return np.vectorize(f)(x)
def vectorize_without_init(x):
return vf(x)
b = perfplot.bench(
setup=np.random.rand,
n_range=[2 ** k for k in range(20)],
kernels=[
f,
array_for,
array_map,
fromiter,
vectorize,
vectorize_without_init,
],
xlabel="len(x)",
)
b.save("out1.svg")
b.show()
使用 numpy.vectorize
怎么样。
import numpy as np
x = np.array([1, 2, 3, 4, 5])
squarer = lambda t: t ** 2
vfunc = np.vectorize(squarer)
vfunc(x)
# Output : array([ 1, 4, 9, 16, 25])
The vectorize function is provided primarily for convenience, not for performance. The implementation is essentially a for loop.
在其他问题中,我发现 vectorize
可能会使用户迭代速度加倍。但真正的加速是使用真正的 numpy
数组操作。
squarer(x)
已经适用于非一维数组。 vectorize
仅比列表理解(如问题中的那个)真正具有任何优势,而不是 squarer(x)
。
np.vectorize
比等效的列表理解要慢。现在它可以更好地扩展,因此使用大参数时它会更快。它仍然不如使用已编译的 numpy
方法和运算符而没有任何类型的 python 级别循环那么快。
TL;博士
正如 @user2357112 所指出的,应用函数的“直接”方法始终是将函数映射到 Numpy 数组的最快和最简单的方法:
import numpy as np
x = np.array([1, 2, 3, 4, 5])
f = lambda x: x ** 2
squares = f(x)
通常避免使用 np.vectorize
,因为它表现不佳,并且有(或有)多个 issues。如果您正在处理其他数据类型,您可能需要研究下面显示的其他方法。
方法比较
这里有一些简单的测试来比较映射函数的三种方法,这个例子使用 Python 3.6 和 NumPy 1.15.4。一、测试的设置功能:
import timeit
import numpy as np
f = lambda x: x ** 2
vf = np.vectorize(f)
def test_array(x, n):
t = timeit.timeit(
'np.array([f(xi) for xi in x])',
'from __main__ import np, x, f', number=n)
print('array: {0:.3f}'.format(t))
def test_fromiter(x, n):
t = timeit.timeit(
'np.fromiter((f(xi) for xi in x), x.dtype, count=len(x))',
'from __main__ import np, x, f', number=n)
print('fromiter: {0:.3f}'.format(t))
def test_direct(x, n):
t = timeit.timeit(
'f(x)',
'from __main__ import x, f', number=n)
print('direct: {0:.3f}'.format(t))
def test_vectorized(x, n):
t = timeit.timeit(
'vf(x)',
'from __main__ import x, vf', number=n)
print('vectorized: {0:.3f}'.format(t))
使用五个元素进行测试(从最快到最慢排序):
x = np.array([1, 2, 3, 4, 5])
n = 100000
test_direct(x, n) # 0.265
test_fromiter(x, n) # 0.479
test_array(x, n) # 0.865
test_vectorized(x, n) # 2.906
拥有 100 多个元素:
x = np.arange(100)
n = 10000
test_direct(x, n) # 0.030
test_array(x, n) # 0.501
test_vectorized(x, n) # 0.670
test_fromiter(x, n) # 0.883
并且有 1000 个或更多的数组元素:
x = np.arange(1000)
n = 1000
test_direct(x, n) # 0.007
test_fromiter(x, n) # 0.479
test_array(x, n) # 0.516
test_vectorized(x, n) # 0.945
不同版本的 Python/NumPy 和编译器优化会有不同的结果,所以对你的环境做一个类似的测试。
count
参数和生成器表达式,那么 np.fromiter
会明显更快。
'np.fromiter((f(xi) for xi in x), x.dtype, count=len(x))'
f(x)
的直接解,which beats everything else by over an order of magnitude。
f
有 2 个变量并且数组是 2D 的呢?
周围有 numexpr、numba 和 cython,这个答案的目标是考虑这些可能性。
但首先让我们明确一点:无论您如何将 Python 函数映射到 numpy 数组,它仍然是一个 Python 函数,这意味着对于每次评估:
numpy-array 元素必须转换为 Python 对象(例如 Float)。
所有计算都使用 Python 对象完成,这意味着需要解释器、动态调度和不可变对象的开销。
因此,由于上面提到的开销,使用哪种机器来实际循环遍历数组并没有起到很大的作用——它比使用 numpy 的内置功能要慢得多。
让我们看一下下面的例子:
# numpy-functionality
def f(x):
return x+2*x*x+4*x*x*x
# python-function as ufunc
import numpy as np
vf=np.vectorize(f)
vf.__name__="vf"
np.vectorize
被选为纯 Python 函数类方法的代表。使用 perfplot
(请参阅此答案附录中的代码),我们得到以下运行时间:
https://i.stack.imgur.com/LZdLl.png
我们可以看到,numpy-approach 比纯 python 版本快 10 到 100 倍。更大数组大小的性能下降可能是因为数据不再适合缓存。
值得一提的是,vectorize
也使用了大量内存,因此内存使用通常是瓶颈(参见相关的 SO-question)。另请注意,numpy 在 np.vectorize
上的文档指出它“主要是为了方便,而不是为了性能”。
应该使用其他工具,当需要性能时,除了从头开始编写 C 扩展外,还有以下可能性:
人们经常听到,numpy 的性能是最好的,因为它在引擎盖下是纯 C。不过还有很大的提升空间!
矢量化的 numpy 版本使用大量额外的内存和内存访问。 Numexp-library 尝试平铺 numpy-arrays,从而获得更好的缓存利用率:
# less cache misses than numpy-functionality
import numexpr as ne
def ne_f(x):
return ne.evaluate("x+2*x*x+4*x*x*x")
导致以下比较:
https://i.stack.imgur.com/uhW6K.png
我无法解释上图中的所有内容:一开始我们可以看到 numexpr-library 的开销更大,但是因为它更好地利用了缓存,所以对于更大的数组,它的速度大约快了 10 倍!
另一种方法是对函数进行 jit 编译,从而获得真正的纯 C UFunc。这是 numba 的方法:
# runtime generated C-function as ufunc
import numba as nb
@nb.vectorize(target="cpu")
def nb_vf(x):
return x+2*x*x+4*x*x*x
它比原来的 numpy-approach 快 10 倍:
https://i.stack.imgur.com/LBPzB.png
然而,这个任务是令人尴尬的并行化,因此我们也可以使用 prange
来并行计算循环:
@nb.njit(parallel=True)
def nb_par_jitf(x):
y=np.empty(x.shape)
for i in nb.prange(len(x)):
y[i]=x[i]+2*x[i]*x[i]+4*x[i]*x[i]*x[i]
return y
正如预期的那样,并行函数对于较小的输入较慢,但对于较大的输入则更快(几乎是因子 2):
https://i.stack.imgur.com/e3djr.png
虽然 numba 专注于使用 numpy-arrays 优化操作,但 Cython 是一个更通用的工具。提取与 numba 相同的性能更为复杂 - 通常取决于 llvm (numba) 与本地编译器 (gcc/MSVC):
%%cython -c=/openmp -a
import numpy as np
import cython
#single core:
@cython.boundscheck(False)
@cython.wraparound(False)
def cy_f(double[::1] x):
y_out=np.empty(len(x))
cdef Py_ssize_t i
cdef double[::1] y=y_out
for i in range(len(x)):
y[i] = x[i]+2*x[i]*x[i]+4*x[i]*x[i]*x[i]
return y_out
#parallel:
from cython.parallel import prange
@cython.boundscheck(False)
@cython.wraparound(False)
def cy_par_f(double[::1] x):
y_out=np.empty(len(x))
cdef double[::1] y=y_out
cdef Py_ssize_t i
cdef Py_ssize_t n = len(x)
for i in prange(n, nogil=True):
y[i] = x[i]+2*x[i]*x[i]+4*x[i]*x[i]*x[i]
return y_out
Cython 会导致一些较慢的功能:
https://i.stack.imgur.com/Z0mEw.png
结论
显然,只测试一个功能并不能证明任何事情。还应该记住,对于选择的函数示例,内存的带宽是大于 10^5 个元素的大小的瓶颈 - 因此我们在该区域对 numba、numexpr 和 cython 具有相同的性能。
最后,最终答案取决于函数类型、硬件、Python 分布和其他因素。例如,Anaconda-distribution 将 Intel 的 VML 用于 numpy 的函数,因此在 exp
、sin
、cos
等先验函数和类似函数中的性能优于 numba(除非它使用 SVML,请参阅此 SO-post) - 参见例如以下SO-post。
然而,根据这次调查和我迄今为止的经验,我想说,只要不涉及超越函数,numba 似乎是最简单且性能最佳的工具。
使用 perfplot-package 绘制运行时间:
import perfplot
perfplot.show(
setup=lambda n: np.random.rand(n),
n_range=[2**k for k in range(0,24)],
kernels=[
f,
vf,
ne_f,
nb_vf, nb_par_jitf,
cy_f, cy_par_f,
],
logx=True,
logy=True,
xlabel='len(x)'
)
squares = squarer(x)
数组上的算术运算是按元素自动应用的,具有高效的 C 级循环,避免了适用于 Python 级循环或理解的所有解释器开销。
您想要应用于 NumPy 数组元素的大多数函数都可以正常工作,尽管有些可能需要更改。例如,if
在元素方面不起作用。您希望将它们转换为使用 numpy.where
之类的结构:
def using_if(x):
if x < 5:
return x
else:
return x**2
变成
def using_where(x):
return numpy.where(x < 5, x, x**2)
似乎没有人提到在 numpy 包中生成 ufunc
的内置工厂方法:np.frompyfunc
,我已经针对 np.vectorize
进行了测试,并且性能优于它约 20~30%。当然,它的性能不如规定的 C 代码甚至 numba
(我没有测试过),但它可以比 np.vectorize
更好的选择
f = lambda x, y: x * y
f_arr = np.frompyfunc(f, 2, 1)
vf = np.vectorize(f)
arr = np.linspace(0, 1, 10000)
%timeit f_arr(arr, arr) # 307ms
%timeit vf(arr, arr) # 450ms
我还测试了更大的样本,并且改进是成比例的。另请参阅文档here
编辑: 原始答案具有误导性, np.sqrt
直接应用于数组,只是开销很小。
在多维情况下,您希望应用对一维数组进行操作的内置函数,numpy.apply_along_axis 是一个不错的选择,也适用于来自 numpy 和 scipy 的更复杂的函数组合。
先前的误导性陈述:
添加方法:
def along_axis(x):
return np.apply_along_axis(f, 0, x)
perfplot 代码的性能结果接近 np.sqrt
。
f
进行矢量化。例如,尝试在 Nico 的 perf 代码中将 np.sqrt
替换为 math.sqrt
,您将收到错误消息。这里实际发生的是使用数组参数调用 f
,因为 x 是一维的,并且您告诉它沿包含所有元素的第一个轴应用它。要使此答案有效,应将 apply_along_axis
的参数替换为 x[None,:]
。然后你会发现,long_axis 是其中最慢的。
np.sqrt
。
我相信在 numpy 的较新版本(我使用 1.13)中,您可以通过将 numpy 数组传递给您为标量类型编写的函数来简单地调用该函数,它会自动将函数调用应用于 numpy 数组上的每个元素并返回给您另一个numpy数组
>>> import numpy as np
>>> squarer = lambda t: t ** 2
>>> x = np.array([1, 2, 3, 4, 5])
>>> squarer(x)
array([ 1, 4, 9, 16, 25])
t
的每个元素 t 的是 **
运算符。那是普通的numpy。将它包装在 lambda
中不会做任何额外的事情。
如 this post 中所述,只需使用生成器表达式,如下所示:
numpy.fromiter((<some_func>(x) for x in <something>),<dtype>,<size of something>)
以上所有答案都比较好,但是如果您需要使用自定义函数进行映射,并且您有numpy.ndarray
,则需要保留数组的形状。
我只比较了两个,但它会保留 ndarray
的形状。我使用了包含 100 万个条目的数组进行比较。这里我使用了 square 函数,它也是 numpy 内置的,并且有很大的性能提升,因为需要一些东西,你可以使用你选择的函数。
import numpy, time
def timeit():
y = numpy.arange(1000000)
now = time.time()
numpy.array([x * x for x in y.reshape(-1)]).reshape(y.shape)
print(time.time() - now)
now = time.time()
numpy.fromiter((x * x for x in y.reshape(-1)), y.dtype).reshape(y.shape)
print(time.time() - now)
now = time.time()
numpy.square(y)
print(time.time() - now)
输出
>>> timeit()
1.162431240081787 # list comprehension and then building numpy array
1.0775556564331055 # from numpy.fromiter
0.002948284149169922 # using inbuilt function
在这里您可以清楚地看到 numpy.fromiter
考虑到简单的方法效果很好,如果内置功能可用,请使用它。
fromiter
快 8% .. 这可能不会改变游戏规则(即可能不值得额外的认知负担)。
使用numpy.fromfunction(function, shape, **kwargs)
请参阅“https://docs.scipy.org/doc/numpy/reference/generated/numpy.fromfunction.html”
function
。这不是OP想要的。
不定期副业成功案例分享
f(x)
排除在您的情节之外。它可能不适用于每个f
,但它适用于此处,并且在适用时它很容易成为最快的解决方案。vf = np.vectorize(f); y = vf(x)
因短输入而获胜的说法。