python - 生成给定总和的随机(int)列表python

Question

我想用给定的概率prob模拟n个随机选择。

我当前的解决方案如下:

from random import choices

result = [0]*len(prob)
population = list(range(0,len(prob)))
ch = choices(population, weights=prob, k=n)
for i in ch:
   result[i] += 1

我的问题是，我多次调用此代码，并且通常使用较大的 n ，并且此解决方案似乎根本没有效率。

是否有更好的方法(例如某些库的预构建函数)？

总而言之，我想要构建一个随机列表的最有效方法，总和为 $n$，这样获得给定列表的概率等于作为 n 个随机选择获得此列表的概率概率prob。

谢谢

[编辑以添加上下文]

我真正做的是在一种马尔可夫链中进行n次随机游走，如下所示:

def rand_walk(n,state):
    next_states,prob = complicated_function(state) // compute the next states and their probability
    succ = distribute_over_next_states(n, prob) // compute how many walk goes to which states
    accu = complicated_function_2(state) // accumulator for the result
    for ns in range(0,len(next_states)):
        accu += rand_walk(succ[i],next_states[I])
    return accu

重点是，下一个状态及其概率的计算成本很高，因此我避免多次计算(因此我避免按顺序运行n)。这就是为什么我想按照给定的概率分配n。

我希望这足以清楚地理解......

Answer 1

哼。 Numpy 已经实现 multinomial draws ，所以我们甚至不需要 via_binomial 函数:

In [56]: np.random.multinomial(10**8, [0.2, 0.3, 0.5])
Out[56]: array([20003098, 29996630, 50000272])

In [57]: via_binomial(10**8, [0.2, 0.3, 0.5])
Out[57]: [19993527, 30000996, 50005477]

---

IIUC，您可以将其视为一系列重复的二项式绘制:

from random import choices
import numpy as np

def original(n, prob):
    result = [0]*len(prob)
    population = list(range(0,len(prob)))
    ch = choices(population, weights=prob, k=n)
    for i in ch:
       result[i] += 1
    return result

def via_binomial(n, prob):
    result = []
    already_handled_prob = 0
    n_left = n
    for p in prob[:-1]:
        draw_prob = p / (1 - already_handled_prob)
        result.append(np.random.binomial(n_left, draw_prob))
        already_handled_prob += p
        n_left -= result[-1]
    result.append(n-sum(result))
    return result

给我

In [29]: %time original(10**6, [1])
Wall time: 343 ms
Out[29]: [1000000]

In [30]: %time via_binomial(10**6, [1])
Wall time: 0 ns
Out[30]: [1000000]

In [31]: %time original(10**6, [0.25, 0.75])
Wall time: 343 ms
Out[31]: [249944, 750056]

In [32]: %time via_binomial(10**6, [0.25, 0.75])
Wall time: 0 ns
Out[32]: [250030, 749970]

In [33]: %time original(10**8, [0.4, 0.3, 0.2, 0.1])
Wall time: 40.1 s
Out[33]: [40004163, 29999878, 19992540, 10003419]

In [34]: %time via_binomial(10**8, [0.4, 0.3, 0.2, 0.1])
Wall time: 0 ns
Out[34]: [39997530, 29999334, 20003182, 9999954]

In [35]: %time via_binomial(10**8, [0.4, 0.3, 0.2, 0.1])
Wall time: 0 ns
Out[35]: [40009324, 29995955, 19996223, 9998498]

(好吧，我使用的是 %time 而不是 %timeit 来作弊。:-) 在我的机器上大约需要 3-10 us。)