python 用于评估随机游走的库?

python library for evaluating random walks?

我正在尝试评估随机游走结束位置的概率,但我在程序速度方面遇到了一些问题。基本上我想做的是将包含随机游走概率的字典作为输入(例如 p = {0:0.5, 1:0.2. -1:0.3} 意思有50%的概率X保持在0,20%的概率X增加1,30%的概率X减少1)然后计算n次迭代后所有可能的未来状态的概率。

例如,如果 p = {0:0.5,1:0.2。 -1:0.3} 和 n = 2 那么它将 return {0:0.37, 1:0.2, -1:0.3, 2:0.04 , -2:0.09} 如果 p = {0:0.5, 1:0.2。 -1:0.3} 和 n = 1 那么它将 return {0:0.5, 1:0.2。 -1:0.3}

我有工作代码,如果 n 很低并且 p 字典很小,它运行得相对较快,但是当 n > 500 并且字典有大约 50 个值时,它需要超过 5 分钟的时间来计算。我猜这是因为它只在一个处理器上执行所以我继续修改它所以它会使用 python 的多处理模块(因为我读到多线程不会提高并行计算性能,因为 GIL ).

我的问题是,多处理没有太大改进,现在我不确定是因为我实施错误还是因为 python 中多处理的开销。我只是想知道是否有某个库可以并行评估 n > 500 时随机游走的所有可能性的所有概率?如果找不到任何东西,我的下一步是在 C 中编写自己的函数作为扩展,但这将是我第一次这样做,尽管我已经用 C 编写了一段时间。

原始非多处理代码

def random_walk_predictor(probabilities_tree, period):
    ret = probabilities_tree
    probabilities_leaves = ret.copy()
    for x in range(period):
        tmp = {}
        for leaf in ret.keys():
            for tree_leaf in probabilities_leaves.keys():
                try:
                    tmp[leaf + tree_leaf] = (ret[leaf] * probabilities_leaves[tree_leaf]) + tmp[leaf + tree_leaf]
                except:
                    tmp[leaf + tree_leaf] = ret[leaf] * probabilities_leaves[tree_leaf]
        ret = tmp
    return ret

多处理代码

from multiprocessing import Manager,Pool
from functools import partial

def probability_calculator(origin, probability, outp, reference):
    for leaf in probability.keys():
        try:
            outp[origin + leaf] = outp[origin + leaf] + (reference[origin] * probability[leaf])
        except KeyError:
            outp[origin + leaf] = reference[origin] * probability[leaf]

def random_walk_predictor(probabilities_leaves, period):
    probabilities_leaves = tree_developer(probabilities_leaves)
    manager = Manager()
    prob_leaves = manager.dict(probabilities_leaves)
    ret = manager.dict({0:1})
    p = Pool()

    for x in range(period):
        out = manager.dict()
        partial_probability_calculator = partial(probability_calculator, probability = prob_leaves, outp = out, reference = ret.copy())

        p.map(partial_probability_calculator, ret.keys())
        ret = out

    return ret.copy()

往往有分析解决方案来精确解决这种看起来类似于二项式分布的 problem,但我假设你真的在寻求更一般的计算解决方案 class 的问题。

与使用 python 词典相比,从潜在的数学问题的角度来考虑这个问题更容易。构建一个矩阵 A 来描述从一种状态转到另一种状态的概率。建立一个状态 x 来描述在某个时间处于给定位置的概率。

因为在 n 转换之后,您最多可以从原点(在任一方向)步进 n 步 - 您的状态需要有 2n+1 行,并且 A 需要为正方形,大小为 2n+1 x 2n+1。

对于两个时间步长的问题,您的转换矩阵将为 5x5,如下所示:

[[ 0.5  0.2  0.   0.   0. ]
 [ 0.3  0.5  0.2  0.   0. ]
 [ 0.   0.3  0.5  0.2  0. ]
 [ 0.   0.   0.3  0.5  0.2]
 [ 0.   0.   0.   0.3  0.5]]

你在时间 0 的状态将是:

[[ 0.]
 [ 0.]
 [ 1.]
 [ 0.]
 [ 0.]]

系统的一步演化可以通过乘以Ax来预测。

所以在 t = 1 时,

 x.T = [[ 0.   0.2  0.5  0.3  0. ]]

并且在 t = 2 时,

x.T = [[ 0.04  0.2   0.37  0.3   0.09]]

因为即使是适度的时间步数,这也可能会占用相当多的存储空间(A 需要 n^2 存储空间),但是非常稀疏,我们可以使用稀疏矩阵来减少我们的存储空间(并加快我们的计算速度)。这样做意味着 A 需要大约 3n 个元素。

import scipy.sparse as sp
import numpy as np

def random_walk_transition_probability(n, left = 0.3, centre = 0.5, right = 0.2):
    m = 2*n+1
    A  = sp.csr_matrix((m, m))
    A += sp.diags(centre*np.ones(m), 0)
    A += sp.diags(left*np.ones(m-1), -1)
    A += sp.diags(right*np.ones(m-1),  1)
    x = np.zeros((m,1))
    x[n] = 1.0
    for i in xrange(n):
        x = A.dot(x)
    return x

print random_walk_transition_probability(4)

计时

%timeit random_walk_transition_probability(500)
100 loops, best of 3: 7.12 ms per loop

%timeit random_walk_transition_probability(10000)
1 loops, best of 3: 1.06 s per loop