4.7 示例:随机漫步

下面是一个通过内置的random模块以纯Python的方式实现1000步的随机漫步:

图4-4是根据前100个随机漫步值生成的折线图:

    1. In [251]: nsteps = 1000
    2. In [252]: draws = np.random.randint(0, 2, size=nsteps)
    3. In [253]: steps = np.where(draws > 0, 1, -1)
    4. In [254]: walk = steps.cumsum()

    有了这些数据之后,我们就可以沿着漫步路径做一些统计工作了,比如求取最大值和最小值:

    现在来看一个复杂点的统计任务——首次穿越时间,即随机漫步过程中第一次到达某个特定值的时间。假设我们想要知道本次随机漫步需要多久才能距离初始0点至少10步远(任一方向均可)。np.abs(walk)>=10可以得到一个布尔型数组,它表示的是距离是否达到或超过10,而我们想要知道的是第一个10或-10的索引。可以用argmax来解决这个问题,它返回的是该布尔型数组第一个最大值的索引(True就是最大值):

    1. In [257]: (np.abs(walk) >= 10).argmax()
    2. Out[257]: 37

    注意,这里使用argmax并不是很高效,因为它无论如何都会对数组进行完全扫描。在本例中,只要发现了一个True,那我们就知道它是个最大值了。

    1. In [258]: nwalks = 5000
    2. In [259]: nsteps = 1000
    3. In [260]: draws = np.random.randint(0, 2, size=(nwalks, nsteps)) # 0 or 1
    4. In [261]: steps = np.where(draws > 0, 1, -1)
    5. In [262]: walks = steps.cumsum(1)
    6. In [263]: walks
    7. Out[263]:
    8. array([[ 1, 0, 1, ..., 8, 7, 8],
    9. [ 1, 0, -1, ..., 34, 33, 32],
    10. [ 1, 0, -1, ..., 4, 5, 4],
    11. ...,
    12. [ 1, 2, 3, ..., 14, 13, 14],
    13. [ -1, -2, -3, ..., -24, -23, -22]])

    现在,我们来计算所有随机漫步过程的最大值和最小值:

    得到这些数据之后,我们来计算30或-30的最小穿越时间。这里稍微复杂些,因为不是5000个过程都到达了30。我们可以用any方法来对此进行检查:

    1. In [266]: hits30 = (np.abs(walks) >= 30).any(1)
    2. In [267]: hits30
    3. Out[267]: array([False, True, False, ..., False, True, False], dtype=bool)
    4. In [268]: hits30.sum() # Number that hit 30 or -30
    5. Out[268]: 3410

    然后我们利用这个布尔型数组选出那些穿越了30(绝对值)的随机漫步(行),并调用argmax在轴1上获取穿越时间:

    1. In [269]: crossing_times = (np.abs(walks[hits30]) >= 30).argmax(1)
    2. In [270]: crossing_times.mean()