2

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

例如，如果 p = {0:0.5, 1:0.2。-1:0.3} 和 n = 2 如果 p = {0:0.5, 1:0.2，它将返回 {0:0.37, 1:0.2, -1:0.3, 2:0.04, -2:0.09}。-1:0.3} 和 n = 1 那么它将返回 {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()
4

1 回答 1

2

往往有分析解决方案来精确解决这种看起来类似于二项分布的问题，但我假设您真的要求为更一般的问题类别提供计算解决方案。

与其使用 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.]]

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

所以在 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
于 2015-07-17T13:33:39.457 回答