在 [0..n-1] 范围内生成 m 个不同的随机数

Posted 2023-02-25

【中文标题】在 [0..n-1] 范围内生成 m 个不同的随机数

【英文标题】：Generating m distinct random numbers in the range [0..n-1]

【发布时间】：2011-10-20 07:53:49

【问题描述】：

我有两种方法可以

方法一：

//C++-ish pseudocode
int result[m];
for(i = 0; i < m; ++i)

   int r;
   do
   
      r = rand()%n;
   while(r is found in result array at indices from 0 to i)
   result[i] = r;   

方法二：

//C++-ish pseudocode
int arr[n];
for(int i = 0; i < n; ++i)
    arr[i] = i;
random_shuffle(arr, arr+n);
result = first m elements in arr;

当 n 远大于 m 时，第一种方法更有效，而第二种方法在其他情况下更有效。但是“更大”并不是一个严格的概念，是吗？ :)

问题： 我应该使用 n 和 m 的什么公式来确定 method1 还是 method2 更有效？ （就运行时间的数学期望而言）

【问题讨论】：

如果m 真的很小，那么效率有那么重要吗？针对更可能导致问题的情况进行优化。

@Mark：我会动态获取n和m。我将不得不确定运行时使用哪种方法

你用几个不同的参数测试过这两种方法吗？只是为了大致了解他们需要多长时间。

【参考方案1】：

纯数学： 让我们计算一下这两种情况下rand()函数调用的数量并比较结果：

案例 1： 当您已经选择了 k 个数字时，让我们看看在步骤 i = k 上调用的数学期望。通过一次rand() 调用获得数字的概率等于p = (n-k)/n。我们需要知道这种调用数量的数学期望，这会导致获得一个我们还没有的数字。

使用1 调用得到它的概率是p。使用 2 调用 - q * p，其中 q = 1 - p。在一般情况下，在n 调用之后恰好得到它的概率是(q^(n-1))*p。因此，数学期望是 Sum[ n * q^(n-1) * p ], n = 1 --> INF。这个总和等于1/p（由 wolfram alpha 证明）。

因此，在步骤i = k 中，您将执行1/p = n/(n-k) 对rand() 函数的调用。

现在让我们总结一下：

Sum[ n/(n - k) ], k = 0 --> m - 1 = n * T - 方法 1 中的rand 调用次数。 这里T = Sum[ 1/(n - k) ], k = 0 --> m - 1

案例 2：

这里rand() 在random_shuffle n - 1 内被调用了多次（在大多数实现中）。

现在，要选择方法，我们必须比较这两个值：n * T ? n - 1。 因此，要选择适当的方法，请按上述计算T。如果T < (n - 1)/n 最好使用第一种方法。否则使用第二种方法。

【讨论】：

如果你能扩展会很棒 - “使用 2 个调用 - q * p，其中 q = 1 - p。”我真的不明白。

我还会考虑 while 循环每次迭代中的搜索工作量（案例 2 为 0），而不仅仅是随机调用的时间

【参考方案2】：

查看original Fisher-Yates algorithm 的***描述。它提倡基本上使用您的方法 1 最多 n/2，其余部分使用您的方法 2。

【讨论】：

我认为这是一个调整问题。他只需要在目标系统上测量它。

这听起来很对，当然从内存效率的角度来看。使用较小的m 值可能会不必要地消耗大量内存。

在visualstudiomagazine.com/articles/2013/07/01/…，您可以找到Fisher 方法的一个很好解释的c# 实现。这篇文章还展示了水库法的另一种方法以及它们之间的区别。

【参考方案3】：

就个人而言，我会使用方法 1，然后如果 M > N/2，则选择 N-M 个值，然后反转数组（返回未选择的数字）。例如，如果 N 为 1000，而您想要其中的 950 个，则使用方法 1 选择 50 个值，然后返回其他 950 个。

编辑：虽然，如果一致的性能是您的目标，我会使用修改后的方法 2，它不会进行完全随机播放，而只会随机播放 N 长度数组的前 M 个元素。

int arr[n];
for(int i = 0; i < n; ++i)
    arr[i] = i;

for (int i =0; i < m; ++i) 
   int j = rand(n-i); // Pick random number from 0 <= r < n-i.  Pick favorite method
   // j == 0 means don't swap, otherwise swap with the element j away
   if (j != 0)  
      std::swap(arr[i], arr[i+j]);
   

result = first m elements in arr;

【参考方案4】：

对于任何结果集，这是一个适用于 O(n) 内存和 O(n) 时间（其中 n 是返回结果的数量，而不是您从中选择的集合的大小）的算法。为了方便，它在 Python 中是因为它使用哈希表：

def random_elements(num_elements, set_size):
    state = 
    for i in range(num_elements):
        # Swap state[i] with a random element
        swap_with = random.randint(i, set_size - 1)
        state[i], state[swap_with] = state.get(swap_with, swap_with), state.get(i, i)
    return [state[i] for i in range(num_elements) # effectively state[:num_elements] if it were a list/array.

这只是一个部分的 Fisher-yates 洗牌，被洗牌的数组被实现为一个稀疏哈希表 - 任何不存在的元素都等于它的索引。我们打乱第一个num_elements 索引，并返回这些值。在set_size = 1, 的情况下，这等效于在范围内选择一个随机数，而在num_elements = set_size 的情况下，这等效于标准的fisher-yates shuffle。

观察到这是 O(n) 时间是微不足道的，因为循环的每次迭代最多初始化哈希表中的两个新索引，所以它也是 O(n) 空间。

【讨论】：

哈希表/字典访问不是常数时间 O(1) 而是对数时间 O(log(n))，因此总体复杂度为 O(n log(n))。

@astraujums 实际上是常数摊销

我认为应该是swap_with = random.randint(i, set_size-1)，因为 randint() 使用包含范围？ @尼克-约翰逊

@cidermole 哇，你是对的。不错的 Python 库，这是我以前从未注意到的。固定。

【参考方案5】：

第三种方法呢？

int result[m];
for(i = 0; i < m; ++i)

   int r;
   r = rand()%(n-i);
   r += (number of items in result <= r)
   result[i] = r;   

---

编辑应该是

这样更好，一个使用来自 Fisher-Yates 的 Modern Method 的示例

//C++-ish pseudocode
int arr[n];
for(int i = 0; i < n; ++i)
    arr[i] = i;

for(i = 0; i < m; ++i)
    swap(arr, n-i, rand()%(n-i) );

result = last m elements in arr;

【讨论】：

(number of items in result < r) 是什么意思？这是否意味着如果 r 大于结果中的项目数，则要添加 1？

+1 再一次，我们错过了显而易见的事情。但它会占用大量空间，因为您还必须将结果保存在排序列表中才能有效地执行(number of items in result < r) 检查。

@Nobody 我的意思是为结果中的每个项目 +1 rand() to n-i

的事实

@Jacob：我的意思是，即使我们假设 rand() 返回一个均匀分布的数字，你的算法也不会生成均匀分布的序列

@Armen 为什么不呢？本质上不是说在[0,n-1] 中选择一个随机整数，然后从剩余的n-1 中选择一个随机整数……最后从剩余的n-m+1 中选择一个。这是Fisher-Yates shuffle 的开始，但在 n-m 而不是 2 处停止。

【参考方案6】：

谈论数学期望，这没什么用，但我还是会发布它：D

Shuffle 是简单的 O(m)。

现在其他算法有点复杂。生成下一个数字所需的步数是试验次数的期望值，试验长度的概率是几何分布。所以...

p=1          E[X1]=1            = 1           = 1
p=1-1/n      E[x2]=1/(1-1/n)    = 1 + 1/(n-1) = 1 + 1/(n-1) 
p=1-2/n      E[x3]=1/(1-1/n)    = 1 + 2/(n-2) = 1 + 1/(n-2) + 1/(n-2)
p=1-3/n      E[X4]=1/(1-2/n)    = 1 + 3/(n-3) = 1 + 1/(n-3) + 1/(n-3) + 1(n-3)
....
p=1-(m-1)/n) E[Xm]=1/(1-(m-1)/n))

请注意，总和可以分成三角形，见右手边。

让我们使用调和级数的公式：H_n = Sum k=0->n (1/k) = approx ln(k)

Sum(E[Xk]) = m + ln(n-1)-ln(n-m-1) + ln(n-2)-ln(n-m-1) + ... = m + ln(n-1) + ln(n-2) + ... - (m-1)*ln(n-m-1) ..

还有一些关于调和级数和的论坛，如果你仍然有兴趣我会查一下...

更新：实际上这是一个相当不错的公式（感谢出色的具体数学书）

Sum(H_k) k=0->n = n*H_n - n

所以预期的步数：

Sum(E[Xk]) = m + (n-1)*ln(n-1) - (n-1) - (n-m-1)*ln(n-m-1) - (n-m-1)) - (m-1)*ln(n-m-1).

注意：我还没有验证。

【参考方案7】：

用set代替数组怎么样，我觉得比数组容易多了

set<int> Numbers;
while (Numbers.size() < m) 
   Numbers.insert(rand() % n);

【参考方案8】：

这有点远，但它可能会起作用，具体取决于您的系统。

从一些合理的比率开始，例如 0.5。 当一个请求进来时，用你从当前阈值比率中得到的任何一种方法来处理它。 记录所花费的时间，当您有“空闲”时间时，使用其他方法执行相同的任务。 如果替代解决方案比原始解决方案快得多，请向上或向下调整阈值。

这种方法的明显缺陷是，在负载高度可变的系统上，您的“离线”测试不会太可靠。

【参考方案9】：

有人建议费舍尔-耶茨洗牌。不知道接下来的代码是否会生成均匀分布的整数，但至少是紧凑且一次性的：

std::random_device rd;
std::mt19937 g(rd());
for (size_type i = 1; i < std::size(v); ++i) 
    v[i] = std::exchange(v[g() % i], i);

【参考方案10】：

很可能在调试模式下启动它（并保留一个方法作为注释）几次以获得平均值会更简单，然后使用另一种方法从中获取平均值

【参考方案11】：

我不建议这种方法，但它有效

#include <iostream>
#include <random>
#include <ctime>

using namespace std;

int randArray[26];
int index = 0;

bool unique(int rand) 

    for (int i = 0; i < index; i++)
        if (rand == randArray[i])
            return false;
    index++;
    return true;



int main()

    srand(time(NULL));

    for (int i = 1; i < 26; i++)
        randArray[i] = -1;

    for (int i = 0; i < 26; i++) 

        randArray[i] = rand() % 26;

        while (!unique(randArray[i])) 
            randArray[i] = rand() % 26;
        
    

    for (int i = 0; i < 26; i++) 
        cout << randArray[i] << " ";
    

    cout << "\n" << index << endl;


    return 0;