2.4. Optimizing code

作者GaëlVaroquaux

本章介绍使Python代码运行速度更快的策略。

先决条件

2.4.1. 优化工作流程

  1. 使其工作:以简单的清晰的方式编写代码。
  2. 使它可靠地工作:编写自动测试用例,确保你的算法是正确的,如果你打破它,测试将捕获破损。
  3. 通过分析简单的用例来找到瓶颈并加快这些瓶颈,找到更好的算法或实现来优化代码。请记住,在一个现实例子的分析和代码的简单性和执行速度之间应该找到一个折衷。为了高效的工作,最好使用大约10秒的分析运行。

2.4.2. 剖析Python代码

无需测量即可优化!

  • 测量:剖析,计时
  • 你会有惊喜:最快的代码并不总是你的想法

2.4.2.1. Timeit

在IPython中,使用timeithttps://docs.python.org/library/timeit.html)来计算基本操作:

In [1]: import numpy as np
In [2]: a = np.arange(1000)
In [3]: %timeit a ** 2
100000 loops, best of 3: 5.73 us per loop
In [4]: %timeit a ** 2.1
1000 loops, best of 3: 154 us per loop
In [5]: %timeit a * a
100000 loops, best of 3: 5.56 us per loop

使用此引导您在策略之间的选择。

注意

对于长时间通话,请使用%time而不是%timeit;它不太精确,但更快

2.4.2.2. Profiler

当您有一个大型程序进行配置文件时,例如following file时,

# For this example to run, you also need the 'ica.py' file
import numpy as np
from scipy import linalg
from ica import fastica
def test():
data = np.random.random((5000, 100))
u, s, v = linalg.svd(data)
pca = np.dot(u[:, :10].T, data)
results = fastica(pca.T, whiten=False)
if __name__ == '__main__':
test()

注意

这是两种无监督学习技术的组合,主成分分析(PCA)和独立成分分析(ICA)。PCA是用于降维的技术,即,使用较少维度来解释数据中观察到的方差的算法。ICA是源分离技术,例如用于解混合通过多个传感器记录的多个信号。如果您有比信号更多的传感器,首先执行PCA,然后执行ICA可能很有用。有关详细信息,请参见:来自scikits-learn的FastICA示例。

要运行它,您还需要下载ica module在IPython中,我们可以对脚本进行计时:

In [1]: %run -t demo.py
IPython CPU timings (estimated):
User : 14.3929 s.
System: 0.256016 s.

并配置它:

In [2]: %run -p demo.py
916 function calls in 14.551 CPU seconds
Ordered by: internal time
ncalls tottime percall cumtime percall filename:lineno (function)
1 14.457 14.457 14.479 14.479 decomp.py:849 (svd)
1 0.054 0.054 0.054 0.054 {method 'random_sample' of 'mtrand.RandomState' objects}
1 0.017 0.017 0.021 0.021 function_base.py:645 (asarray_chkfinite)
54 0.011 0.000 0.011 0.000 {numpy.core._dotblas.dot}
2 0.005 0.002 0.005 0.002 {method 'any' of 'numpy.ndarray' objects}
6 0.001 0.000 0.001 0.000 ica.py:195 (gprime)
6 0.001 0.000 0.001 0.000 ica.py:192 (g)
14 0.001 0.000 0.001 0.000 {numpy.linalg.lapack_lite.dsyevd}
19 0.001 0.000 0.001 0.000 twodim_base.py:204 (diag)
1 0.001 0.001 0.008 0.008 ica.py:69 (_ica_par)
1 0.001 0.001 14.551 14.551 {execfile}
107 0.000 0.000 0.001 0.000 defmatrix.py:239 (__array_finalize__)
7 0.000 0.000 0.004 0.001 ica.py:58 (_sym_decorrelation)
7 0.000 0.000 0.002 0.000 linalg.py:841 (eigh)
172 0.000 0.000 0.000 0.000 {isinstance}
1 0.000 0.000 14.551 14.551 demo.py:1 (<module>)
29 0.000 0.000 0.000 0.000 numeric.py:180 (asarray)
35 0.000 0.000 0.000 0.000 defmatrix.py:193 (__new__)
35 0.000 0.000 0.001 0.000 defmatrix.py:43 (asmatrix)
21 0.000 0.000 0.001 0.000 defmatrix.py:287 (__mul__)
41 0.000 0.000 0.000 0.000 {numpy.core.multiarray.zeros}
28 0.000 0.000 0.000 0.000 {method 'transpose' of 'numpy.ndarray' objects}
1 0.000 0.000 0.008 0.008 ica.py:97 (fastica)
...

显然,svd(在decomp.py)是什么需要我们的大部分时间,a.k.a.瓶颈。我们必须找到一种方法使这一步更快,或避免这一步(算法优化)。花费时间对剩余的代码是没用的。

在IPython外部进行概要分析,运行“cProfile”

可以在IPython外部进行类似的分析,只需调用内置的Python分析器 cProfileprofile

$  python -m cProfile -o demo.prof demo.py

使用-o开关将会将分析器结果输出到文件demo.prof,以使用外部工具查看。如果您希望使用可视化工具处理分析器输出,这将非常有用。

2.4.2.3. 线轮廓线

分析器告诉我们哪个函数占用大部分时间,但不是在哪里被调用。

为此,我们在源文件中使用line_profiler:来装饰我们要用@profile检查的几个函数(不需要导入它)

@profile
def test():
data = np.random.random((5000, 100))
u, s, v = linalg.svd(data)
pca = np.dot(u[:, :10], data)
results = fastica(pca.T, whiten=False)

然后使用kernprof.py程序运行脚本,使用开关-l, - 逐行 t1>和-v, - view使用逐行分析器并查看结果,

$ kernprof.py -l -v demo.py
Wrote profile results to demo.py.lprof
Timer unit: 1e-06 s
File: demo.py
Function: test at line 5
Total time: 14.2793 s
Line # Hits Time Per Hit % Time Line Contents
=========== ============ ===== ========= ======= ==== ========
5 @profile
6 def test():
7 1 19015 19015.0 0.1 data = np.random.random((5000, 100))
8 1 14242163 14242163.0 99.7 u, s, v = linalg.svd(data)
9 1 10282 10282.0 0.1 pca = np.dot(u[:10, :], data)
10 1 7799 7799.0 0.1 results = fastica(pca.T, whiten=False)

SVD正在占用所有时间。我们需要优化这一线。

2.4.3. 使代码更快速

一旦我们确定了瓶颈,我们需要使相应的代码更快。

2.4.3.1. 算法优化

首先要寻找的是算法优化:有没有方法来计算更少,或更好?

对于问题的高级视图,对算法背后的数学的良好理解有帮助。然而,常见的是发现简单的变化,例如在for循环之外移动计算或内存分配,这带来了巨大的收益。

2.4.3.1.1. SVD示例

在上面的两个例子中,SVD - 奇异值分解 - 是大部分时间所需要的。实际上,该算法的计算成本在输入矩阵的大小中大致为n^3

但是,在这两个示例中,我们不使用SVD的所有输出,而是仅使用其第一个返回参数的前几行。如果我们使用scipy的svd实现,我们可以要求SVD的不完整版本。注意,线性代数在scipy中的实现比那些在numpy中更丰富,应该是首选。

In [3]: %timeit np.linalg.svd(data)
1 loops, best of 3: 14.5 s per loop
In [4]: from scipy import linalg
In [5]: %timeit linalg.svd(data)
1 loops, best of 3: 14.2 s per loop
In [6]: %timeit linalg.svd(data, full_matrices=False)
1 loops, best of 3: 295 ms per loop
In [7]: %timeit np.linalg.svd(data, full_matrices=False)
1 loops, best of 3: 293 ms per loop

然后,我们可以使用此洞察optimize the previous code

def test():
data = np.random.random((5000, 100))
u, s, v = linalg.svd(data, full_matrices=False)
pca = np.dot(u[:, :10].T, data)
results = fastica(pca.T, whiten=False)
In [1]: import demo
In [2]: %timeit demo.
demo.fastica demo.np demo.prof.pdf demo.py demo.pyc
demo.linalg demo.prof demo.prof.png demo.py.lprof demo.test
In [2]: %timeit demo.test()
ica.py:65: RuntimeWarning: invalid value encountered in sqrt
W = (u * np.diag(1.0/np.sqrt(s)) * u.T) * W # W = (W * W.T) ^{-1/2} * W
1 loops, best of 3: 17.5 s per loop
In [3]: import demo_opt
In [4]: %timeit demo_opt.test()
1 loops, best of 3: 208 ms per loop

实际不完全SVD,例如。只计算前10个特征向量,可以使用arpack计算,可在scipy.sparse.linalg.eigsh中获得。

计算线性代数

对于某些算法,许多瓶颈将是线性代数计算。在这种情况下,使用正确的函数来解决正确的问题是关键。例如,与一般矩阵相比,对称矩阵的特征值问题更容易求解。此外,最常见的情况是,您可以避免反转矩阵,并使用成本更低(且数值更稳定)的操作。

知道你的计算线性代数。遇到疑问时,请浏览scipy.linalg,并使用%timeit尝试不同的替代方法。

2.4.4. 写入更快的数字代码

关于numpy的高级使用的完整讨论在Advanced NumPy章节中或者在van der Walt等人的文章The NumPy array:a structure for efficient numerical computation中找到。这里我们只讨论一些常见的技巧,使代码更快。

  • 向量化for循环

    找到避免使用numpy数组的循环的技巧。为此,掩码和索引数组可能是有用的。

  • 广播

    使用broadcasting在组合数组之前对数组执行尽可能小的操作。

  • 原位操作

    In [1]: a = np.zeros(1e7)
    
    In [2]: %timeit global a ; a = 0*a
    10 loops, best of 3: 111 ms per loop
    In [3]: %timeit global a ; a *= 0
    10 loops, best of 3: 48.4 ms per loop

    注意:我们需要在timeit中使用全局a,以使其工作,因为它分配到a,因此将其视为局部变量。

  • 内存容易:使用视图,而不是副本

    复制大数组与在其上进行简单的数值操作一样昂贵;

    In [1]: a = np.zeros(1e7)
    
    In [2]: %timeit a.copy()
    10 loops, best of 3: 124 ms per loop
    In [3]: %timeit a + 1
    10 loops, best of 3: 112 ms per loop
  • 小心缓存效果

    当分组时,内存访问更便宜:以连续方式访问大数组比随机访问快得多。这暗示了更小的步幅更快(参见CPU cache effects):

    In [1]: c = np.zeros((1e4, 1e4), order='C')
    
    In [2]: %timeit c.sum(axis=0)
    1 loops, best of 3: 3.89 s per loop
    In [3]: %timeit c.sum(axis=1)
    1 loops, best of 3: 188 ms per loop
    In [4]: c.strides
    Out[4]: (80000, 8)

    这就是为什么Fortran排序或C排序可能在操作上产生巨大差异的原因:

    In [5]: a = np.random.rand(20, 2**18)
    
    In [6]: b = np.random.rand(20, 2**18)
    In [7]: %timeit np.dot(b, a.T)
    1 loops, best of 3: 194 ms per loop
    In [8]: c = np.ascontiguousarray(a.T)
    In [9]: %timeit np.dot(b, c)
    10 loops, best of 3: 84.2 ms per loop

    请注意,复制数据以解决此问题可能不值得:

    In [10]: %timeit c = np.ascontiguousarray(a.T)
    
    10 loops, best of 3: 106 ms per loop

    使用numexpr可用于自动优化此类效果的代码。

  • 使用编译代码

    最后一个手段,一旦你确定所有的高级优化已经被探索,是将热点,即几个行或函数,其中大部分时间花费到编译代码。对于编译代码,首选选项是使用Cython:在编译代码中很容易转换退出的Python代码,并且很好地使用numpy支持生成高效代码对numpy数组,例如通过展开循环。

警告

对于所有以上:配置文件和时间你的选择。不要将优化基于理论考虑。