如何有效地在区间内生成排序均匀分布的值?

机器算法验证 随机生成
2022-03-28 22:15:09

假设我想从区间生成一组随机数(a, b)生成的序列还应该具有已排序的属性。我可以想到两种方法来实现这一点。

n为要生成的序列的长度。

第一种算法:

Let `offset = floor((b - a) / n)`
for i = 1 up to n:
   generate a random number r_i from (a, a+offset)
   a = a + offset
   add r_i to the sequence r

第二种算法:

for i = 1 up to n:
    generate a random number s_i from (a, b)
    add s_i to the sequence s
sort(r)

我的问题是,算法 1 生成的序列是否与算法 2 生成的序列一样好?

3个回答

第一个算法严重失败有两个原因:

  1. 的地板可以大大减少它。实际上,当时,它将为零,给您一个值都相同的集合!(ab)/nba<n

  2. 当您不发言时,结果值分布过于均匀 例如,在 iid 均匀变量的任何简单随机样本中(比如在之间),有的机会最大不会在从的上区间。使用算法 1,最大值有的机会出现在该区间内。出于某些目的,这种超均匀性是好的,但总的来说,这是一个可怕的错误,因为 (a) 许多统计数据会被破坏,但 (b) 很难确定原因。na=0b=1(11/n)n1/e37%11/n1100%

  3. 如果您想避免排序,请改为生成独立的指数分布变量。通过除以总和将它们的累积总和标准化为范围删除最大值(始终为)。重新缩放到范围n+1(0,1)1(a,b)

显示了所有三种算法的直方图。独立集合的累积结果,每个值。)算法 1 的直方图中没有任何可见的变化表明存在问题。其他两种算法的变化正是可以预期的——以及随机数生成器所需要的。1000n=100

有关模拟独立均匀变量的更多(有趣)方法,请参阅使用正态分布的绘图模拟均匀分布的绘图

图:直方图

这是R生成该图的代码。

b <- 1
a <- 0
n <- 100
n.iter <- 1e3

offset <- (b-a)/n
as <- seq(a, by=offset, length.out=n)
sim.1 <- matrix(runif(n.iter*n, as, as+offset), nrow=n)
sim.2 <- apply(matrix(runif(n.iter*n, a, b), nrow=n), 2, sort)
sim.3 <- apply(matrix(rexp(n.iter*(n+1)), nrow=n+1), 2, function(x) {
  a + (b-a) * cumsum(x)[-(n+1)] / sum(x)
})

par(mfrow=c(1,3))
hist(sim.1, main="Algorithm 1")
hist(sim.2, main="Algorithm 2")
hist(sim.3, main="Exponential")

第一个算法产生的数字太均匀

另请参阅低差异系列

中的 2 个随机数对于真实的统一数据,机会是 50:50,它们同时大于或小于 0.5。用你的方法,机会是0。所以你的数据统一。[0;1]

(正如所指出的,这可能是分层所需的属性。像 Halton 和 Sobel 这样的低差异系列确实有它们的用例。)

一种适当但昂贵的方法(对于实际价值)

... 是使用 beta 分布的随机数。均匀分布的排序统计量是 beta 分布的。您可以使用它来随机绘制最小的,然后是第二小的,...重复。

假设要在中生成数据。最小值为分布。(对于后续情况,减少并重新调整到剩余间隔)。要生成一般的 beta 随机数,我们需要生成两个 Gamma 分布的随机值。但是然后为此,我们可以从这个分布中抽取随机数作为[0;1]Beta[1,n]n1XBeta[n,1]ln(1X)Exponential[n]ln(U[0;1])n

ln(1x)=ln(1u)n1x=u1nx=1u1n

产生以下算法:

x = a
for i in range(n, 0, -1):
    x += (b-x) * (1 - pow(rand(), 1. / i))
    result.append(x) 

可能涉及数值不稳定性,并且pow每个对象的计算和除法可能会比排序慢。

对于整数值,您可能需要使用不同的分布。

排序非常便宜,所以只需使用它

但不要打扰。排序是如此的便宜,所以只需排序。多年来,我们已经很好地理解了如何实现排序双精度数不值得避免的排序算法。理论上它是,但在一个好的实现中,常数项是如此之小,以至于这是一个完美的例子,即理论复杂性结果是多么无用运行基准测试。生成 100 万个有排序和没有排序的随机数。运行几次,如果排序经常超过非排序,我不会感到惊讶,因为排序的成本仍然会比你的测量误差小得多。O(nlogn)

这还取决于您对随机数所做的事情。对于数值积分问题,方法一(当通过去除地板算子进行校正时)将产生优越的点集。你正在做的是一种分层抽样的形式,它的优点是可以避免聚集。例如,不可能在 0-(ba)/n 范围内获得所有值。也就是说,对于其他应用程序,这可能非常糟糕,这取决于你想用它做什么。