Python高速化を試してみた 前半：単純なコードで比較する

※本稿は，[ビジネスデータサイエンス研究会](https://statsforcom.connpass.com)開催の[『ハイパフォーマンスPython 』勉強会【#3】](https://statsforcom.connpass.com/event/110359/)における発表内容の一部を整理し記事化したものです。

# 内容
本稿では，オライリーから出ている[ハイパフォーマンスPython](https://www.amazon.co.jp/dp/4873117402/ref=cm_sw_r_tw_dp_U_x_XrjiCbWE1X5SY )の`7章 Cにコンパイルする`, `9章 multiprocessingモジュール`の内容をもとに，CythonとCPUレベルでの並列処理に絞っていくつかのPython高速化手法を試してみた結果をまとめます。

まず前半では，実用性は皆無だけどとりあえず簡単なコードを題材にして比較を行います。
後半では記事を改めてもう少し実用性のあるコードの高速化を行います。

ソースコードは[こちら](
https://github.com/tatamiya/HPPython_Chapter7_9/tree/master/simple_sum)にアップしております。

# 実行環境

**作業環境（Dockerコンテナ）**
OS : Ubuntu 18.04.1
CPU数 : 4
メモリ容量 : 16GB

**ホストマシン**
OS : Windows10 Home Edition
CPU数 : 8
メモリ容量 : 16GB

# 対象アルゴリズム
1~Nまでの和を，N=1, 2, 3, ..., 100000について別々に計算したリストを作成します。

```Python:modules/function.py
def sum_from_one(N):
    sum = 0
    for i in range(1, N+1):
        sum += i

    return sum
```

```Python:simple_sum_normal.py
import numpy as np
from modules import function

if __name__ == "__main__":
    sum_list = []
    repeat_list = np.random.permutation(np.arange(1, 100001))
    for i in repeat_list:
        sum_list.append(function.sum_from_one(i))
```

実行時間は1200s程度でした。

### 比較内容

以下の実装で実行速度，改変のしやすさの比較を行いました。

- そのままのPython
- Cython
- 並列処理（CPUレベル）
    - multiprocessing
    - joblib
- Cython + 並列処理（OpenMP）

[元の本](https://www.amazon.co.jp/dp/4873117402/ref=cm_sw_r_tw_dp_U_x_XrjiCbWE1X5SY )では他にも`Numba`の利用や`CPython`実装の改変，プロセス間通信込みの並列処理なども紹介されていましたが，今回はなるべく既存のコードへの修正が少なくて済む範囲ということで，上記に絞って調べました。
他の高速化手法や，本ではあえて触れられていなかったGPUによる並列化についても，機会があれば調査して記事を執筆したいと考えております。

# Cython

単純に1~Nを足し合わせる関数`sum_from_one(N)`の部分だけ変数の型を定義し，CythonでC言語にコンパイルして読み込みました。

```Python:cython_function.pyx
def sum_from_one_cython(int N):
    cdef unsigned int sum
    cdef unsigned int i

    sum = 0

    for i in range(1, N+1):
        sum += i

    return sum
```


```Python:simple_sum_cython.py
import numpy as np
from modules import calculate

if __name__ == "__main__":
    sum_list = []
    repeat_list = np.random.permutation(np.arange(1, 100001))
    for i in repeat_list:
        sum_list.append(calculate.sum_from_one_cython(i))
```

### コンパイル方法

Cythonはコードをそのまま実行することができず，一旦コンパイルしてからライブラリとして読み込む必要があります。

そのために，コンパイル用スクリプト`setup.py`を作成してからこちらを実行しました。

```Python:setup.py
from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext

setup(
    cmdclass={'build_ext': build_ext},
    ext_modules=[Extension("calculate", ["cython_function.pyx"]
)
```

```
# python setup.py build_ext --inplace
```

なお，読み込み用のライブラリは`calculate`と名付けてあります。
実行本体の`.py`ファイルでは，こちらを`import calculate`として読み込み，関数`calculate.sum_from_one_cython()`を呼び出して計算を行います。

### 実行結果

```
In [1] %run -t simple_sum_cython.py
IPython CPU timings (estimated):
  User   :       0.48 s.
  System :       0.00 s.
Wall time:       0.48 s.
```

型を定義しただけですが，素のPythonのときより2000倍以上速くなっています。

# 並列処理

[ハイパフォーマンスPython](https://www.amazon.co.jp/dp/4873117402/ref=cm_sw_r_tw_dp_U_x_XrjiCbWE1X5SY ) 9章で紹介されているmultiprocessingモジュールだけでなく，9.4節末尾のコラムで触れられているjoblibも使ってみました。

### **multiprocessing**

```Python:simple_sum_multiprocess.py
from multiprocessing import Pool
import multiprocessing as multi
import numpy as np
from modules import function

if __name__ == "__main__":

    repeat_list = np.random.permutation(np.arange(1, 100001))

    p = Pool(multi.cpu_count())
    p.map(function.sum_from_one, repeat_list)
    p.close()
```

```:実行結果
In [1]: %run -t simple_sum_multiprocess.py

IPython CPU timings (estimated):
  User   :       0.99 s.
  System :       0.06 s.
Wall time:     299.23 s.
```

### **joblib**

```Python:simple_sum_joblib.py
from joblib import Parallel, delayed
import numpy as np
from modules import function

if __name__ == "__main__":

    repeat_list = np.random.permutation(np.arange(1, 100001))

    r = Parallel(n_jobs=-1)([delayed(function.sum_from_one)(i) for i in repeat_list])
```

```:実行結果
In [3]: %run -t simple_sum_joblib.py

IPython CPU timings (estimated):
  User   :      14.53 s.
  System :       0.38 s.
Wall time:     296.49 s.
```

### 並列処理考察

素のPythonでの実行時間が1400秒だったことを思い出すと，ちょうどコア数4個分の1/4程度に短縮しています。ただし，Cythonの足元にはまだ及びません。
`multiprocessing`と`joblib`ではほとんど実行時間が変わりませんが，コードの行数で見ると`joblib`のほうが少なく，内容も若干ですが理解しやすいように思えます。

ちなみに`joblib`はscikit-learnの実装でも使われています。
（例）[GridSearchCV](https://github.com/scikit-learn/scikit-learn/blob/master/sklearn/model_selection/_search.py)
扱いやすさが評価されているのかもしれません。

# Cython + OpenMP

Cythonを使うことで圧倒的に速くなり，また並列処理もコア数の分だけ高速化されていたので，この２つを組み合わせればさらなる高速化が望めるだろうということでOpenMPを使った高速化も試してみました。

### コード

今回のコードについて言えば，for文を回すところで`range`を`prange`に置き換えるだけです。

```Python:cython_function.pyx
#cython: boundscheck=False
from cython.parallel import prange

def sum_from_one(int N):
    cdef unsigned int sum
    cdef unsigned int i

    sum = 0
    for i in prange(1, N+1, nogil=True):
        sum += i

    return sum
```

なお，`prange`の中に`nogil=True`というオプションがあります。
これは，**G**lobal **I**nterpreter **L**ock，通称GILを解除する，という意味です。
通常PythonではGILによって単一スレッドでしかプログラムが実行されないように制限されていますが，これを外すことで複数コア・スレッドでの処理を陽に指定しています。

### コンパイル

基本的には先ほどの単純Cythonと一緒ですが，`setup.py`に追記して，`-fopenmp`を指定する必要があります。
詳しくは[Gitにあげたコード](https://github.com/tatamiya/HPPython_Chapter7_9/blob/master/simple_sum/modules/setup.py)をご参照ください。

### 実行結果

```Python
import numpy as np
import calculate

if __name__ == "__main__":
    sum_list = []
    repeat_list = np.random.permutation(np.arange(1, 100001))
    for i in repeat_list:
        sum_list.append(calculate.sum_from_one(i))
```

```:実行結果
In [1]: %run -t simple_sum_openmp.py

IPython CPU timings (estimated):
  User   :       5.38 s.
  System :       0.93 s.
Wall time:       1.61 s.
```

素のPythonや`joblib`, `multiprocessing`による並列処理と比べると圧倒的に早いのですが，単純なCythonからは3倍遅くなってしまっています。

### 小さな処理は並列化すると逆に遅くなる

並列化しているのになんで遅くなるの？と疑問に思うかもしれませんが，上記では足しあげる数のMaxを$N=10^5$でやっていましたが，これを変えてみると以下のようになります。

|繰り返し回数$N$|単純Cython|OpenMP|
|---|---|---|
|$10^4$|0.02s|0.13s|
|$10^5$|0.48s|1.61s|
|$10^6$|43.95s|46.69s|

$N$が大きくなるにつれて，型指定だけのCythonとOpenMPを使った場合との差が小さくなっていくのがわかります。
並列処理を行う場合，*複数コアに分散させる処理のかたまりがある程度大きくないと，どこかしらでオーバーヘッドがかかって遅くなってしまう*，ということが考えられます。
なお，同様の現象は`joblib`, `multiprocessing`による並列化でも見られました。

# まとめと考察

以上の内容を表にまとめると，以下の通りです。

|実装内容|実行時間|扱いやすさ|
|---|---|---|
|Python|1178s|普通に書くだけ|
|Cython（型定義）|0.48s|ビルドの手間がかかる|
|並列化(multiprocessing)|299s|少し書き換える|
|並列化(joblib)|296s|少し書き換える|
|Cython + 並列化(OpenMP)|1.61s|ビルドの手間+大きな処理向け|

#### Cythonは速いが手間が増える
ある程度簡単な処理であれば，Cythonで変数の型を指定してやるだけでも圧倒的に速くなるようです。
ただし，ビルドの手間も増えてしまうので，開発のスピードは落ちます。
また，今回触れませんでしたが環境構築も手間がかかる点にも要注意です。特にWindowsでの導入については，少々凶悪だったので[別で記事を作成ました。](https://qiita.com/tatamiya/items/510752d0d7eac5b71ceb)

#### 並列処理は個々のプロセスの大きさに注意
並列処理については，Cythonには及ばないものの`joblib`を使えば比較的楽に導入し高速化することができます。
ただし，あまりに小さな処理だと逆に遅くなるので，並列させる各処理がある程度まとまった大きさを持っている必要があります。これについては，その都度速度比較したうえで選択しましょう。

# おわりに

単純なプログラムを例に高速化手法を試してみましたが，このくらい小さくて簡単な処理であれば，Cythonで型指定するのが最も効率がいいように思えます。

しかし大抵はCythonを導入しなくても，アルゴリズムを改良するか，numpyなど既存パッケージを使った方が，手っ取り早いし余計なビルドの手間もかからないかと思います。

実際今回の例についても，そもそもN-1での結果にNを足して更新していくか，もしくは式
$$ \frac{N(N+1)}{2}$$
を使って内包表記で書けば済む話でした。

```Python
In [1]: [N * (N + 1) / 2 for N in range(1,100001)]
```

次回は，もう少し複雑で実用性もある問題について高速化を試していきたいと思います。