我想在 python 中测试向量化代码的性能:

import timeit
import numpy as np

def func1():
  x = np.arange(1000)
  sum = np.sum(x*2)
  return sum

def func2():
  sum = 0
  for i in xrange(1000):
    sum += i*2
  return sum

def func3():
  sum = 0
  for i in xrange(0,1000,4):
    x = np.arange(i,i+4,1)
    sum += np.sum(x*2)
  return sum

print timeit.timeit(func1, number = 1000)
print timeit.timeit(func2, number = 1000)
print timeit.timeit(func3, number = 1000)

代码给出以下输出:

0.0105729103088
0.069864988327
0.983253955841

第一个和第二个功能的性能差异并不奇怪。但令我惊讶的是，第三个功能明显比其他功能慢。

与 Python 相比，我对 C 语言的矢量化代码更加熟悉，并且第三个函数更像 C - 运行一个 for 循环并在每个循环的一条指令中处理 4 个数字。据我了解，numpy 调用 C 函数，然后将 C 中的代码矢量化。因此，如果是这种情况，我的代码也一次将 4 个数字传递给 numpy。当我一次传递更多数字时，代码的性能不会更好。那么为什么它要慢得多呢？是因为调用 numpy 函数的开销吗？

此外，我什至首先想出第三个函数的原因是因为我担心 func1 中 x 的大量内存分配的性能。

我的担心有道理吗？为什么以及如何改进它或为什么不改进它？

提前致谢。

编辑:

出于好奇，虽然它违背了我创建第 3 个版本的初衷，但我研究了 roganjosh 的建议并尝试了以下编辑。

def func3():
  sum = 0
  x = np.arange(0,1000)
  for i in xrange(0,1000,4):
    sum += np.sum(x[i:i+4]*2)
  return sum

输出:

0.0104308128357
0.0630609989166
0.748773813248

有提高，但与其他功能相比仍有较大差距。

是不是因为x[i:i+4]还是新建了一个数组？

编辑 2:

我已经按照大牛的建议再次修改了代码。

def func1():
  x = np.arange(1000)
  x *= 2
  return x.sum()

def func3():
  sum = 0
  x = np.arange(0,1000)
  for i in xrange(0,1000,4):
    x[i:i+4] *= 2
    sum += x[i:i+4].sum()
  return sum

输出:

0.00824999809265
0.0660569667816
0.598328828812

还有另一个加速。所以numpy数组的声明肯定是个问题。现在func3中应该只有一个数组声明，但时间还是慢了很多。是因为调用 numpy 数组的开销吗？

最佳答案

您似乎最感兴趣的是函数 3 与纯 NumPy(函数 1)和 Python(函数 2)方法相比的区别。答案很简单(尤其是如果您查看函数 4):

• NumPy 函数有一个“巨大的”常数因子。

您通常需要数千个元素才能进入 np.sum 的运行时间实际上取决于数组中元素数量的状态。使用 IPython 和 matplotlib(图在答案末尾)，您可以轻松检查运行时依赖性:

import numpy as np

n = []
timing_sum1 = []
timing_sum2 = []
for i in range(1, 25):
    num = 2**i
    arr = np.arange(num)
    print(num)
    time1 = %timeit -o arr.sum()    # calling the method
    time2 = %timeit -o np.sum(arr)  # calling the function
    n.append(num)
    timing_sum1.append(time1)
    timing_sum2.append(time2)

np.sum(缩写)的结果非常有趣:

4
22.6 µs ± 297 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)
16
25.1 µs ± 1.08 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
64
25.3 µs ± 1.58 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
256
24.1 µs ± 1.48 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
1024
24.6 µs ± 221 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)
4096
27.6 µs ± 147 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)
16384
40.6 µs ± 1.29 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
65536
91.2 µs ± 1.03 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
262144
394 µs ± 8.09 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
1048576
1.24 ms ± 4.38 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
4194304
4.71 ms ± 22.9 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
16777216
18.6 ms ± 280 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

在我的计算机上，常数因子似乎大约是 20µs)并且它需要一个包含 16384 千个元素的数组才能使该时间翻倍。所以函数 3 和 4 的时序主要是常数因子的时序乘法。

在函数 3 中，您包含了 2 次常数因子，一次使用 np.sum，一次使用 np.arange。在这种情况下，arange 非常便宜，因为每个数组的大小都相同，因此 NumPy 和 Python 以及您的操作系统可能会重用最后一次迭代的数组的内存。然而，即使这样也需要时间(对于我计算机上非常小的阵列，大约 2µs)。

更一般地说:要确定瓶颈，您应该始终分析功能!

我用 line-profiler 显示函数的结果.因此我稍微改变了函数，让它们每行只做一个操作:

import numpy as np

def func1():
    x = np.arange(1000)
    x = x*2
    return np.sum(x)

def func2():
    sum_ = 0
    for i in range(1000):
        tmp = i*2
        sum_ += tmp
    return sum_

def func3():
    sum_ = 0
    for i in range(0, 1000, 4):  # I'm using python3, so "range" is like "xrange"!
        x = np.arange(i, i + 4, 1)
        x = x * 2
        tmp = np.sum(x)
        sum_ += tmp
    return sum_

def func4():
    sum_ = 0
    x = np.arange(1000)
    for i in range(0, 1000, 4):
        y = x[i:i + 4]
        y = y * 2
        tmp = np.sum(y)
        sum_ += tmp
    return sum_

结果:

%load_ext line_profiler

%lprun -f func1 func1()
Line #      Hits         Time  Per Hit   % Time  Line Contents
==============================================================
     4                                           def func1():
     5         1           62     62.0     23.8      x = np.arange(1000)
     6         1           65     65.0     24.9      x = x*2
     7         1          134    134.0     51.3      return np.sum(x)

%lprun -f func2 func2()
Line #      Hits         Time  Per Hit   % Time  Line Contents
==============================================================
     9                                           def func2():
    10         1            7      7.0      0.1      sum_ = 0
    11      1001         2523      2.5     30.9      for i in range(1000):
    12      1000         2819      2.8     34.5          tmp = i*2
    13      1000         2819      2.8     34.5          sum_ += tmp
    14         1            3      3.0      0.0      return sum_

%lprun -f func3 func3()
Line #      Hits         Time  Per Hit   % Time  Line Contents
==============================================================
    16                                           def func3():
    17         1            7      7.0      0.0      sum_ = 0
    18       251          909      3.6      2.9      for i in range(0, 1000, 4):
    19       250         6527     26.1     21.2          x = np.arange(i, i + 4, 1)
    20       250         5615     22.5     18.2          x = x * 2
    21       250        16053     64.2     52.1          tmp = np.sum(x)
    22       250         1720      6.9      5.6          sum_ += tmp
    23         1            3      3.0      0.0      return sum_

%lprun -f func4 func4()
Line #      Hits         Time  Per Hit   % Time  Line Contents
==============================================================
    25                                           def func4():
    26         1            7      7.0      0.0      sum_ = 0
    27         1           49     49.0      0.2      x = np.arange(1000)
    28       251          892      3.6      3.4      for i in range(0, 1000, 4):
    29       250         2177      8.7      8.3          y = x[i:i + 4]
    30       250         5431     21.7     20.7          y = y * 2
    31       250        15990     64.0     60.9          tmp = np.sum(y)
    32       250         1686      6.7      6.4          sum_ += tmp
    33         1            3      3.0      0.0      return sum_

我不会详细介绍结果，但如您所见，np.sum 绝对是 func3 和 func4。在我写答案之前，我已经猜到 np.sum 是瓶颈，但这些行分析实际验证它是瓶颈。

这导致使用 NumPy 时一个非常重要的事实:

• 知道何时使用它!小型阵列不值得(大部分)。
• 了解 NumPy 函数并直接使用它们。他们已经使用(如果可用的话)编译器优化标志来展开循环。

如果您真的认为某些部分太慢，那么您可以使用:

• NumPy 的 C API 并使用 C 处理数组(使用 Cython 非常容易，但您也可以手动完成)
• Numba(基于 LLVM)。

但一般来说，对于中等大小(几千个条目或更多)的数组，您可能无法击败 NumPy。

时间的可视化:

%matplotlib notebook

import matplotlib.pyplot as plt

# Average time per sum-call
fig = plt.figure(1)
ax = plt.subplot(111)
ax.plot(n, [time.average for time in timing_sum1], label='arr.sum()', c='red')
ax.plot(n, [time.average for time in timing_sum2], label='np.sum(arr)', c='blue')
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_xlabel('elements')
ax.set_ylabel('time it takes to sum them [seconds]')
ax.grid(which='both')
ax.legend()

# Average time per element
fig = plt.figure(1)
ax = plt.subplot(111)
ax.plot(n, [time.average / num for num, time in zip(n, timing_sum1)], label='arr.sum()', c='red')
ax.plot(n, [time.average / num for num, time in zip(n, timing_sum2)], label='np.sum(arr)', c='blue')
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_xlabel('elements')
ax.set_ylabel('time per element [seconds / element]')
ax.grid(which='both')
ax.legend()

这些图是双对数图，我认为这是可视化数据的最佳方式，因为它扩展了几个数量级(我只是希望它仍然可以理解)。

第一个图显示了执行sum需要多少时间:

第二张图显示了执行sum 所需的平均时间除以数组中的元素数。这只是解释数据的另一种方式:

关于python - numpy 中不同矢量化方法的性能，我们在Stack Overflow上找到一个类似的问题： https://stackoverflow.com/questions/43823473/