如何排除 runif() 的特定结果？答案

【问题标题】：How to exclude specific result of runif()?如何排除 runif() 的特定结果？
【发布时间】：2018-06-30 05:05:44
【问题描述】：

当我想在排除特定值（例如 0.5）的情况下在特定间隔内生成带有 runif() 的随机数时，我可以编写这个函数 ex.runif() 谁来完成这项工作，但它比它慢数百倍正常的runif()。谁能指出更好的解决方案？

ex.runif <- function(n, excl, min, max) {
  # ex.runif() excludes the specific value 'excl'
  q <- excl
  while (q == excl) {
    q <- runif(n, min = min, max = max)
  } 
  return(q)
}

set.seed(42)
ex.runif(1, .5, .25, .75)  # exclude .5, interval [.25, .75]
# [1] 0.707403

library(microbenchmark)
microbenchmark(ex.runif(1, .5, .25, .75), runif(1, min = .25, max = .75))
# Unit: microseconds
# expr         min      lq      mean   median       uq     max neval cld
# ex.runif 692.439 704.685 721.51135 715.2735 722.9275 962.373   100   b
# runif      2.041   2.551   3.49044   2.8070   3.3170  21.176   100  a

【问题讨论】：

循环后返回q
这样做的动机是什么？在概率论中，很少有理由区分两个随机变量，它们的支持度相差一组测度为零。
@JohnColeman 我正在从模拟数据中以概率 pi 计算掩码估计量。 pi 被另一个具有概率 p1 的值覆盖（在模拟中随机绘制）。因此 p1 != excl，否则这是不可能的。
Draw random numbers from the mixture of uniform distributions的可能重复

标签： r function random

【解决方案1】：

如果您要排除的值集是有限的，那么在大多数情况下，不需要这样的函数。原因是均匀分布是连续的，并且任何有限数量的值都以概率为零。也就是说，根据概率论，q == excl 的概率为零。

例如，

set.seed(42)
ex.runif(5, .5, .25, .75)
# [1] 0.7074030 0.7185377 0.3930698 0.6652238 0.5708728
set.seed(42)
runif(5, 0.25, 0.75)
# [1] 0.7074030 0.7185377 0.3930698 0.6652238 0.5708728

同样的情况也很可能在任何其他种子下发生。因此，您可以继续使用runif。

@duckmayr 很好地说明了数值精度。事实上，随着[min, max] 的间隔越来越窄，q == excl 成为真的概率越来越高，在某些应用中，它甚至可能变得相关。

但是，如果理论上您确实只需要排除单个值 0.5，那么执行像 q == excl 这样的检查甚至可能会因为排除不必要的抽奖而造成伤害。

例如，在我的情况下，.Machine$double.eps 是 2.220446e-16。那么当[min,max] 是[0.5 - 10^(-k), 0.5 + 10^(-k)] 并且得出错误结论时，从[0.5 - .Machine$double.eps / 4, 0.5 + .Machine$double.eps / 4] 获得平局的概率是 2 * (2.220446e-16 / 4) / (2 * 10^(-k)) 或大约 0.55 * 10^(k-16)。

【讨论】：

但是，在 0.5 左右的（任意小的）区间内绘制数字的概率为正；例如，在 [0.5 - .Machine$double.eps/4, 0.5 + .Machine$double.eps/4] 内（特别重要，因为 0.5 == 0.5 + .Machine$double.eps / 4 返回 TRUE）。当然，这个概率本身非常小（在我的机器上大约为 0.00000000000000006）。
punif(0.5 + .Machine$double.eps/4) - punif(0.5 - .Machine$double.eps/4)，我会通过sprintf()：sprintf("%0.17f", punif(0.5 + .Machine$double.eps/4) - punif(0.5 - .Machine$double.eps/4))
一个快速实验来强调遇到 0.5 的可能性有多大：x <- runif(10**8) 之后评估 which(x==0.5) 没有发现任何命中。
限制碰巧只是有害的，如果有的话，因为虽然q == excl 可能成为真的，但它会因为错误的原因（概率为 1）。
@Julius -- 好点子。当然，这意味着对于 OP，如果他们在q == excl 返回TRUE 的情况下进行平局，那么 OP实际上在他们的数据中具有与排除等效的值，所以你' d 有一个问题，您是否错误地丢弃了一个仅因为数值精度而错误地似乎等同于限制的值，或者如果您正确地保留了它，但从您可以测量的情况来看，您的数据中会有一个不可能的值。这似乎是一个非常不可能的问题，但如果它弹出来，就不可能令人满意地解决。