使用 numpy 矢量化或映射来加速循环 - Python NumPy 3D 矩阵 "get rid of a loop" Python 问题，蒙特卡洛

Question

现在我有一个循环填充3D NumPy矩阵。我不太擅长理解3D数组结构，尽管我知道它实际上只是我习惯在2D中思考的正常XxY的XxYxZ表示。如果你想知道这是什么它是一个布朗桥(BB)结构用于蒙特卡罗模拟金融问题。原始代码的功劳(来自作者Kenta Oono在这里修复原始帖子的评论):https://gist.github.com/delta2323/6bb572d9473f3b523e6e。你真的不需要知道它背后的数学原理;它基本上只是分割了一个步骤的路径(在这个例子中是21)，从0开始，应用正态分布的冲击(因此是np.random.randn)，直到它到达最后，也就是0。每条路径都应用于一个模拟价格，随机"冲击"它。随着时间的推移，生成资产在到期时可能遵循的路径。虽然这些是完全不相关的，所以我想我也会传递一个V矩阵来关联路径是正确的，但是，让我们保持简单:

import numpy as np
from matplotlib import pyplot
import timeit
steps = 21
underlyings = 3
sims = 131072
seed = 0 # fix the seed for replicating results
np.random.seed(seed)
def sample_path_batches(underlyings, steps, sims):
dt = 1.0 / (steps-1)
dt_sqrt = np.sqrt(dt)
B = np.empty((underlyings, steps, sims), dtype=float)
B[:,0, :] = 0 # set first step to 0
for n in range(steps - 2):
t = n * dt
xi = np.random.randn(underlyings, sims) * dt_sqrt
B[:, n + 1, :] = B[:, n, :] * (1 - dt / (1 - t)) + xi
B[:, -1, :] = 0 # set last step to 0
return B
start_time = timeit.default_timer()
B = sample_path_batches(underlyings, steps, sims)
print('n' + 'Run time for ', sims, ' simulation steps * underlyings: ', 
np.round((timeit.default_timer() - start_time),3), ' seconds')
pyplot.plot(B[:,:,np.random.randint(0,sims)].T); # plot a random simulation set of paths
pyplot.show()

运行时间为131072个模拟步骤*底层:2.014秒

无论如何，这对我的应用程序来说太慢了，尽管我的第二次内循环的原始版本大约是15秒。我已经看到有人把NumPy矢量化到np。矢量化或使用的地图以"平坦化";一个循环，但我不知道怎么做。我正在寻找一个最佳的"原生python";实现将产生相同的数字。B是3D NumPy数组。如果你想的话，你可以复制粘贴并在线运行:https://mybinder.org/v2/gh/jupyterlab/jupyterlab-demo/HEAD?urlpath=lab/tree/demo

任何建议都是感激的!!即使它只是像这样重构循环，然后应用np.vectorize"不管怎样，我很擅长接受建议，并从一个简单的"新视角"出发，使其发挥作用。如何将问题形象化。我通常只会在Cython (nogil/OpenMP/prange)中做这种事情，但我想知道如何"扁平化"。一般来说是一个循环，在NumPy或Pandas中内置了正常的数学库。

Answer 1

加速这段代码的一个简单解决方案是将并行化它使用Numba。您只需要为功能sample_path_batches使用装饰器@nb.njit('float64[:,:,::1](int64, int64, int64)', parallel=True)(其中nb是Numba模块)。注意，函数中的dtype=float必须替换为dtype=np.float64，这样Numba才能正确编译代码。注意，parallel=True应该自动并行化np.random.randn调用以及循环中的基本后续操作。在10核机器上，快了7倍(使用Numpy需要0.253秒，并行实现Numba需要0.036秒)。如果您没有看到任何改进，您也可以尝试使用prange手动并行化它。

另外，您可以使用np.float32类型可显著提高性能(理论上可提高2倍)。然而，Numpy目前不支持np.random.randn的这种类型。相反，应该使用np.random.default_rng().random(size=underlyings*sims, dtype=np.float32).reshape(underlyings, sims)。不幸的是，Numba可能还不支持它，因为Numpy最近才添加了这个…

如果你有NvidiaGPU另一种解决方案是使用CUDA在GPU上执行该功能。这应该快得多。请注意，Numba有特定的优化功能，可以在使用CUDA的GPU上生成随机的np.float32值(见这里)。