使用 SciPy 或 NumPy 生成具有指定权重的离散随机变量答案

【问题标题】：Generating Discrete random variables with specified weights using SciPy or NumPy使用 SciPy 或 NumPy 生成具有指定权重的离散随机变量
【发布时间】：2012-07-07 13:28:07
【问题描述】：

我正在寻找一个简单的函数，它可以根据它们对应的（也指定的）概率生成指定随机值的数组。我只需要它来生成浮点值，但我不明白为什么它不能生成任何标量。我可以想出许多从现有函数构建它的方法，但我想我可能只是错过了一个明显的 SciPy 或 NumPy 函数。

例如：

>>> values = [1.1, 2.2, 3.3]
>>> probabilities = [0.2, 0.5, 0.3]
>>> print some_function(values, probabilities, size=10)
(2.2, 1.1, 3.3, 3.3, 2.2, 2.2, 1.1, 2.2, 3.3, 2.2)

注意：我找到了 scipy.stats.rv_discrete，但我不明白它是如何工作的。具体来说，我不明白这（下）是什么意思，也不应该做什么：

numargs = generic.numargs
[ <shape(s)> ] = ['Replace with resonable value', ]*numargs

如果我应该使用 rv_discrete，您能否提供一个简单的示例以及对上述“形状”语句的解释？

【问题讨论】：

标签： python random numpy scipy

【解决方案1】：

从离散分布中绘图直接内置于 numpy. 该函数称为random.choice（如果不参考 numpy 文档中的离散分布，很难找到）。

elements = [1.1, 2.2, 3.3]
probabilities = [0.2, 0.5, 0.3]
np.random.choice(elements, 10, p=probabilities)

【讨论】：

太棒了！但是，正确的语法是：np.random.choice(elements, 10, p=list(probabilities))
不错。我认为这个版本是在我发布我的原始问题后发布的（我认为这是在 1.7.0 中首次发布的，我相信它是在 2013 年发布的）。
非常好！似乎也可以在不强制列表的情况下工作：np.random.choice(elements, 10, p=probabilities))。
除了Sina 和zeycus 的cmets 之外，elements 和probabilites 本来可以是普通的lists 而不是numpy.arrays，并且代码的工作方式相同。

【解决方案2】：

这是一个简短的、相对简单的返回加权值的函数，它使用 NumPy 的 digitize、accumulate 和 random_sample。

import numpy as np
from numpy.random import random_sample

def weighted_values(values, probabilities, size):
    bins = np.add.accumulate(probabilities)
    return values[np.digitize(random_sample(size), bins)]

values = np.array([1.1, 2.2, 3.3])
probabilities = np.array([0.2, 0.5, 0.3])

print weighted_values(values, probabilities, 10)
#Sample output:
[ 2.2  2.2  1.1  2.2  2.2  3.3  3.3  2.2  3.3  3.3]

它是这样工作的：

首先使用accumulate 创建垃圾箱。
然后我们使用random_sample创建一堆随机数（在0和1之间）
我们使用digitize 来查看这些数字属于哪些垃圾箱。
并返回相应的值。

【讨论】：

是的，这基本上就是我的想法，但我只是认为可能有一个内置函数可以做到这一点。从它的声音来看，没有这样的事情。我必须承认——我不会做得那么优雅。 - 谢谢
NumPy直接提供numpy.cumsum()，可以用np.add.accumulate()代替（np.add()不是很常用，所以我推荐使用cumsum()）。
+1 为有用的numpy.digitize()！然而，SciPy 实际上提供了一个直接回答问题的函数——见我的回答。
PS：... 正如 Tim_Y 所指出的，使用 SciPy 的函数比使用“手动”解决方案（在 10k 元素上）要慢得多。
这个概率需要标准化吗？

【解决方案3】：

您正朝着一个好的方向前进：内置的scipy.stats.rv_discrete() 非常直接地创建了一个离散随机变量。以下是它的工作原理：

>>> from scipy.stats import rv_discrete  

>>> values = numpy.array([1.1, 2.2, 3.3])
>>> probabilities = [0.2, 0.5, 0.3]

>>> distrib = rv_discrete(values=(range(len(values)), probabilities))  # This defines a Scipy probability distribution

>>> distrib.rvs(size=10)  # 10 samples from range(len(values))
array([1, 2, 0, 2, 2, 0, 2, 1, 0, 2])

>>> values[_]  # Conversion to specific discrete values (the fact that values is a NumPy array is used for the indexing)
[2.2, 3.3, 1.1, 3.3, 3.3, 1.1, 3.3, 2.2, 1.1, 3.3]

上面的分布distrib 因此从values 列表中返回索引。

更一般地，rv_discrete() 在其values=(…,…) 参数的第一个元素中采用一系列整数值，并在这种情况下返回这些值；无需转换为特定的（浮点）值。这是一个例子：

>>> values = [10, 20, 30]
>>> probabilities = [0.2, 0.5, 0.3]
>>> distrib = rv_discrete(values=(values, probabilities))
>>> distrib.rvs(size=10)
array([20, 20, 20, 20, 20, 20, 20, 30, 20, 20])

其中（整数）输入值以所需概率直接返回。

【讨论】：

注意：我尝试在它上面运行 timeit，它似乎比 fraxel 的纯 numpy 版本慢了 100 倍。你知道这是为什么吗？
哇，有趣！在 10k 个元素上，我什至会慢 300 倍。我快速浏览了代码：执行了许多检查，但我想它们无法解释运行时间上的如此大的差异；我对 Scipy 代码的研究不够深入，无法看到差异可能来自哪里……
@TimY 我的天真猜测是，速度慢是因为在纯 Python 中完成了更多的工作，而在 C 中完成了更少的工作（在引擎盖下）。（Python 中的数学/科学包倾向于包装 C 代码。）
假设我从一个概率分布方程开始。不得不用它来为每个值生成一个概率，将其提供给rv_discrete，然后从rv_discrete 返回我开始时分布的近似值，这似乎很愚蠢。有什么方法可以直接使用scipy 使用用户定义的方程？
@dbliss 现在我看到您想到了具有 infinite 个可能值的离散分布的情况（不适合这个问题）。 rv_discrete() 没有这个选项。我不确定这样做的标准方法是什么。（我只能想到将均匀随机变量转换为具有非均匀分布的变量的常用方法的稍微复杂的变体，其中累积概率仅针对最常见的值计算，并在需要时扩展。）

【解决方案4】：

您还可以使用Lea，这是一个专门用于离散概率分布的纯 Python 包。

>>> distrib = Lea.fromValFreqs((1.1,2),(2.2,5),(3.3,3))
>>> distrib
1.1 : 2/10
2.2 : 5/10
3.3 : 3/10
>>> distrib.random(10)
(2.2, 2.2, 1.1, 2.2, 2.2, 2.2, 1.1, 3.3, 1.1, 3.3)

等等！

【讨论】：

【解决方案5】：

最简单的 DIY 方法是将概率总结为累积分布。这样，您将单位间隔拆分为长度等于原始概率的子间隔。现在在 [0,1) 上生成一个统一的随机数，并查看它落在哪个区间。

【讨论】：

欣赏这种数学方式，减少对 python 包的依赖。