从自定义分布生成随机数

Question

我正在尝试从自定义分布中生成随机数，我已经发现了这个问题： Simulate from an (arbitrary) continuous probability distribution 但不幸的是，它对我没有帮助，因为那里建议的方法需要一个分布函数的公式。我的分布是多个均匀分布的组合，基本上分布函数看起来像一个直方图。一个例子是：

f(x) = { 
    0     for  x < 1
    0.5   for  1 <= x < 2
    0.25  for  2 <= x < 4
    0     for  4 <= x
}

Answer 1

你只需要逆CDF方法：

samplef <- function (n) {
  x <- runif(n)
  ifelse(x < 0.5, 2 * x + 1, 4 * x)
  }

自己计算 CDF 以验证：

F(x) = 0                 x < 1
       0.5 * x - 0.5     1 < x < 2
       0.25 * x          2 < x < 4
       1                 x > 4

所以它的倒数是：

invF(x) = 2 * x + 1      0 < x < 0.5
          4 * x          0.5 < x < 1

Answer 2

您可以结合使用从discrete distributions 采样的各种有效方法和连续均匀。

也就是说，从变量的整数部分 Y=[X] 进行模拟，该变量具有离散分布，其概率等于处于每个区间的概率（例如通过表格方法 - 也称为别名方法），然后简单地添加一个随机统一 [0,1$, X = Y+U.

在您的示例中，您让 Y 以 0.5,0.25 和 0.25 的概率取值 1,2,3（这相当于以相等的概率对 1,1,2,3 进行抽样），然后添加一个随机统一。

如果您的“直方图”非常大，这可能是一种非常快速的方法。

在 R 中，您可以通过

做一个简单（如果不是特别有效）的版本

sample(c(1,1,2,3))+runif(1)

或

sample(c(1,1,2,3),n,replace=TRUE)+runif(n)

更一般地，您可以使用sample 中的概率权重参数。

如果您需要比这更快的速度（并且对于某些应用程序，您可能需要，特别是对于大直方图和非常大的样本量），您可以使用链接中提到的方法来加快离散部分的速度，并编程该函数的主力部分，使用较低级别的语言（例如 C 语言）。

也就是说，即使只是使用上面的代码和一个相当“大”的直方图——几十到几百个箱——这种方法似乎——即使在我相当不起眼的笔记本电脑上——也能够很好地生成一百万个随机值不到一秒钟，所以对于许多应用程序来说，这会很好。