この記事は、インテル® デベロッパー・ゾーンに公開されている「Faster random number generation in Intel® Distribution for Python*」(https://software.intel.com/en-us/blogs/2016/06/15/faster-random-number-generation-in-intel-distribution-for-python) の日本語参考訳です。
インテル® Distribution for Python* 2017 Beta Update 1 では、Numpy* の拡張として、numpy.random の設計を反映し、パフォーマンスを向上するためインテル® MKL のベクトル統計ライブラリーを使用した numpy.random_intel が追加されました。
numpy.random を numpy.random_intel に置換するだけで、パフォーマンスを向上できます。
In [1]: import numpy as np In [2]: from numpy import random, random_intel In [3]: %timeit np.random.rand(10**5) 1000 loops, best of 3: 1.05 ms per loop In [4]: %timeit np.random_intel.rand(10**5) 10000 loops, best of 3: 146 µs per loop
出力および時間測定に使用したシステム構成:
ソフトウェア: インテル® Distribution for Python* 3.5.1 2017 Beta Update 1 (Jun. 8, 2016)、Numpy* 1.11.0、インテル® MKL 11.3.3。
ハードウェア: インテル® Xeon® プロセッサー E5-2698 v3 @ 2.30GHz (2 ソケット、16 コア、ハイパースレッディング無効)、64GB RAM。
オペレーティング・システム: Ubuntu* 14.04 LTS。
次の表は、一般的な分布におけるパフォーマンスの向上を示しています。
| 分布 | timing(random) (秒) | timing(random_intel) (秒) | スピードアップ (倍) |
|---|---|---|---|
| uniform(-1, 1) | 0.357 | 0.034 | 10.52 |
| normal(0, 1) | 0.834 | 0.081 | 10.35 |
| gamma(5.2, 1) | 1.399 | 0.267 | 5.25 |
| beta(0.7, 2.5) | 3.677 | 0.556 | 6.61 |
| randint(0, 100) | 0.228 | 0.053 | 4.33 |
| poisson(7.6) | 2.990 | 0.052 | 57.44 |
| hypergeometric(214, 97, 83) | 11.353 | 0.517 | 21.96 |
numpy.random_intel モジュールを使用すると、インテル® MKL でサポートされている、(RandomState クラス・コンストラクターの brng 引数、または seed 初期化メソッドを使用して指定できる) 異なる基本擬似乱数ジェネレーターを利用することができます。デフォルトの基本擬似乱数ジェネレーターは、numpy.random と同じです。
次のコードは、異なる基本擬似乱数ジェネレーターで 100,000 の標準一様確率変数のサンプリングのパフォーマンスを測定し、メルセンヌツイスター基本乱数ジェネレーター ‘MT19937’ と比較します。
import numpy as np
from timeit import timeit, Timer
from numpy import random_intel
from operator import itemgetter
def timer_brng(brng_name):
return Timer('numpy.random_intel.uniform(0,1,size=10**5)',
setup='import numpy.random_intel; numpy.random_intel.seed(77777, brng="{}")'.format(brng_name))
_brngs = ['WH', 'PHILOX4X32X10', 'MT2203', 'MCG59', 'MCG31', 'MT19937', 'MRG32K3A', 'SFMT19937', 'R250']
tdata = sorted([(brng, timer_brng(brng).timeit(number=1000)) for brng in _brngs], key=itemgetter(1))
def relative_perf(tdata):
base = dict(tdata).get('MT19937')
return [(name, t/base) for name, t in tdata]
relative_perf(tdata)
出力は次のとおりです。
In [9]: relative_perf(tdata)
Out[9]:
[('MCG31', 0.607988362606216),
('SFMT19937', 0.6520599397085953),
('MCG59', 0.6522106463148241),
('MT2203', 0.8041550837721791),
('MT19937', 1.0),
('PHILOX4X32X10', 1.021985386332785),
('R250', 1.2991300341128584),
('WH', 1.4285221807604567),
('MRG32K3A', 2.148446327408357)]
numpy.random_intel は、leapfrog および skipahead メソッドを使用して RandomState クラスを初期化し、異なるスレッドで使用されるランダムストリームが相関していないことを保証します。
次に例を示します。
rs1 = np.random_intel.RandomState(77777, brng='SFMT19937') rs2 = np.random_intel.RandomState(77777, brng='SFMT19937') # skip the first stream 2**60 steps ahead rs1.skipahead(2**60) s1 = rs1.randint(0,2**31,size=10**5) s2 = rs2.randint(0,2**31,size=10**5)
ピアソンの相関係数を使用して、2 つのサンプルが相関していないという帰無仮説をテストします。
from scipy.stats.stats import pearsonr pearsonr(s1, s2)
出力は次のとおりです。
In[12]: pearsonr(s1, s2) Out[12]: (0.0034069712002612325, 0.28131564831676692)
データが仮説と矛盾しないことを示しています。
詳細についてはモジュールのドキュメントを、ベクトル統計についてはインテル® MKL のドキュメントをそれぞれ参照してください。
皆様からのフィードバックをお待ちしております。
コンパイラーの最適化に関する詳細は、最適化に関する注意事項を参照してください。

