最有效的方法映射函数在numpy数组

在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])

然而,这可能非常低效,因为我在将新数组转换回numpy数组之前,使用列表推导式将其构造为Python列表。我们能做得更好吗?

914961 次浏览

正如在这篇文章中提到的,只需像这样使用生成器表达式:

numpy.fromiter((<some_func>(x) for x in <something>),<dtype>,<size of something>)

使用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])
squares = squarer(x)

数组上的算术运算以元素方式自动应用,高效的c级循环避免了适用于python级循环或理解的所有解释器开销。

您希望应用到NumPy数组elementwise的大多数函数都可以正常工作,尽管有些函数可能需要更改。例如,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)

博士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,因为它的性能不佳,并且已经(或曾经)使用过许多问题。如果您正在处理其他数据类型,您可能需要研究下面所示的其他方法。

方法比较

下面是一些简单的测试,比较三种映射函数的方法,本例使用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和编译器优化会有不同的结果,所以对您的环境进行类似的测试。

我相信在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])

我已经测试了所有建议的方法,加上np.array(list(map(f, x)))perfplot(我的一个小项目)。

消息#1:如果可以使用numpy的本机函数,就使用它。

如果你试图向量化的函数已经向量化了(就像原始帖子中的x**2例子),使用它是比其他任何东西都快(注意对数尺度):

enter image description here

如果你真的需要向量化,用哪个变量并不重要。

enter image description here


代码重现图:

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()

numexprnumbacython,这个答案的目的是考虑这些可能性。

但首先让我们声明一个显而易见的事实:无论你如何将一个Python函数映射到numpy-array上,它始终是一个Python函数,这意味着对于每一个求值:

  • numpy-array元素必须转换为Python-object(例如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(参见答案附录中的代码),我们得到以下运行时间:

enter image description here

我们可以看到,numpy-方法比纯python版本快10 -100倍。较大数组的性能下降可能是因为数据不再适合缓存。

值得一提的是,vectorize也使用了大量内存,因此内存使用通常是瓶颈(参见相关的问题)。还要注意,numpy在np.vectorize上的文档声明它“主要是为了方便,而不是为了性能”。

如果需要性能,除了从头开始编写c扩展,还可以使用其他工具:


人们经常听说,它的性能非常好,因为它的底层是纯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")

导致以下比较:

enter image description here

我不能在上面的图中解释所有的事情:我们可以看到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方法快10倍:

enter image description here


然而,这个任务是尴尬的可并行的,因此我们也可以使用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倍):

enter image description here


numba专门用于优化numpy-arrays的操作,而Cython是一个更通用的工具。提取与numba相同的性能更加复杂——通常是llvm (numba) vs本地编译器(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会导致一些较慢的函数:

enter image description here


结论

显然,只测试一个函数并不能证明任何东西。还应该记住,对于所选的函数(例如,内存带宽对于大于10^5个元素的大小是瓶颈),因此我们在这个区域中对numba、numexpr和cython具有相同的性能。

最后,最终的答案取决于函数类型、硬件、python分布和其他因素。例如,Anaconda-distribution对numpy的函数使用英特尔的VML,因此对于超越函数,如expsincos和类似的函数,它的性能很容易超过numba(除非它使用SVML,请参阅位置),例如下面的位置

然而,从这次调查和到目前为止我的经验来看,我想说,只要不涉及先验函数,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)'
)

似乎没有人提到在numpy包中生成ufunc的内置工厂方法:np.frompyfunc,我已经对np.vectorize进行了测试,并且性能优于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

我也测试了更大的样本,改进是成比例的。请参见文档在这里

# eyz1 # eyz2 # eyz0 # eyz3。

在多维情况下,您希望应用一个操作1d数组的内置函数,numpy.apply_along_axis是一个不错的选择,对于numpy和scipy的更复杂的函数组合也是如此。

先前的误导性陈述:

添加方法:

def along_axis(x):
return np.apply_along_axis(f, 0, x)

perfplot代码给出接近np.sqrt的性能结果。

使用# EYZ0

看到“# EYZ0”

以上答案比较好,但如果你需要使用自定义函数进行映射,你有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工作得很好,考虑到简单的方法,如果内置函数可用,请使用它。