0

我想使用 R 随机生成一个整数序列，每个整数都是从整数池（0,1,2,3....,k）中挑选出来的，并进行替换。k 是预先确定的。(0,1,2,3....,k) 中每个整数 k 的选择概率是 p k (1-p)，其中 p 是预先确定的。也就是说，与 k 相比，选择 1 的概率要高得多，并且我的最终整数序列的 1 可能比 k 多。我不确定如何在 R 中实现这个数字选择过程。

4

2 回答 2

2

此类问题的通用方法是：

  1. 计算p^k * (1-p)每个整数的
  2. 在表中创建这些的累积总和t。
  3. 从均匀分布中画一个数字range(t)
  4. 测量t该数字下降的距离并检查对应的整数。
  5. 整数的概率越大，它将覆盖的范围越大。

这是快速而肮脏的示例代码：

draw <- function(n=1, k, p) {
    v <- seq( 0, k )
    pr <- (p ** v) * (1-p)
    t <- cumsum(pr)
    r <- range(t)
    x <- runif( n, min=min(r), max=max(r) )
    f <- findInterval( x, vec=t )
    v[ f+1 ] ## first interval is 0, and it will likely never pass highest interval
}

请注意，建议的解决方案并不关心您的密度函数加起来是否为 1。根据您的描述，在现实生活中它可能会。但这对于解决方案并不重要。

于 2021-03-13T19:06:54.973 回答
1

天狼星的回答很好。但正如我所知，您所描述的类似于截断的几何分布。

我应该注意到几何分布在不同的作品中定义不同（例如参见MathWorld），所以我们使用如下定义的分布：

  • P(X = x) ~ p^x * (1 - p)，其中 x 是 [0, k] 中的整数。

我对 R 不是很熟悉，但解决方案涉及调用rgeom(1, 1 - p)直到结果为k或更少。

或者，您可以使用通用拒绝采样器，因为概率是已知的（这里更好地称为权重，因为它们不需要总和为 1）。拒绝采样描述如下：

将权重存储在列表中。计算最高权重，称其为max。k然后，使用拒绝采样在区间 [0, ] 中选择一个整数：

  1. i在区间 [0, k]中选择一个均匀的随机整数。
  2. 以概率weights[i]/max（weights[i] = p^i * (1-p)在你的情况下），返回i。否则，转到步骤 1。

给定每个项目的权重，除了拒绝抽样或 Sirius 答案中的解决方案外，还有许多其他方法可以做出加权选择；请参阅我关于加权选择算法的说明。

于 2021-03-13T21:14:05.873 回答