本文用Python模拟随机漫步行为。

1 使用内建的的random模块

import random

position = 0
walk = [position]
steps = 1000
for i in range(steps):
    step = 1 if random.randint(0, 1) else -1
    position += step
    walk.append(position)

random模块每次只能生成一个样本值,效率很低。如果要生成大量样本值,可用numpy.random模块。

可用下面的代码测试两者生成\(1,000,000\)个样本值的速度:

import numpy as np
from random import normalvariate

N = 1000000

%timeit samples = [normalvariate(0, 1) for _ in range(N)]
%timeit np.random.normal(size=N)

输出:

1.17 s ± 14.6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
34.1 ms ± 504 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

2 使用numpy.random

import numpy as np

np.random.seed(12345)

nsteps = 1000
draws = np.random.randint(0, 2, size=nsteps)
steps = np.where(draws > 0, 1, -1)
walk = steps.cumsum()

注意,random.randint(a,b)函数生成的是\([a,b]\)之间的整数,而numpy.random.randint(a,b)生成的是\([a,b)\)之间的整数。

最终生成的walknumpy.ndarray类型的数据。可以看一下最小和最大到达的地方,以及在何时首先偏离原点10的距离。

print(walk.min())
print(walk.max())
print((np.abs(walk) >= 10).argmax())

输出:

-3
31
37

3 同时模拟多个随机漫步

使用numpy.random()可以同时模拟多个随机漫步。这里同时模拟\(5,000\)个,步长依旧设为\(1,000\)

nwalks = 5000
nsteps = 1000
draws = np.random.randint(0, 2, size=(nwalks, nsteps)) # 0 or 1
steps = np.where(draws > 0, 1, -1)
walks = steps.cumsum(axis=1)

得到的walks,是一个\(5,000\times 1,000\)的矩阵,也是numpy.ndarray类型的数据。

一共有多少次随机漫步,达到过偏离原点30的距离?

hits30 = (np.abs(walks) >= 30).any(1)
hits30
hits30.sum() # Number that hit 30 or -30

输出:

3412

在这些随机漫步过程中,平均用了多少步才偏离原点30?

crossing_times = (np.abs(walks[hits30]) >= 30).argmax(1)
crossing_times.mean()

输出:

497.04103165298943