本文用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)\)之间的整数。
最终生成的walk
为numpy.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