生成 numpy.random.choice 的二维数组，无需替换

Question

我很想通过删除一些 for 循环和使用数组来加快我的代码速度。现在最慢的一步是随机列表的生成。

上下文：我的染色体中有许多突变，我想执行 1000 条具有相同长度和相同突变数量的随机“染色体”，但它们的位置是随机的。

这是我目前正在运行以生成这些随机突变位置：

iterations=1000
Chr_size=1000000
num_mut=500
randbps=[]
for k in range(iterations):
listed=np.random.choice(range(Chr_size),num_mut,replace=False)
randbps.append(listed)

我想做一些类似于他们在this question 中介绍的事情

np.random.choice(range(Chr_size),size=(num_mut,iterations),replace=False)

但不替换则适用于整个数组。

进一步的上下文：稍后在脚本中，我遍历每个随机染色体并计算给定窗口中的突变数量：

for l in range(len(randbps)):

    arr=np.asarray(randbps[l])

    for i in range(chr_last_window[f])[::step]:

        counter=((i < arr) & (arr < i+window)).sum()

Answer 1

我不知道 np.random.choice 是如何实现的，但我猜它已针对一般情况进行了优化。另一方面，您的数字不太可能产生相同的序列。在这种情况下，从头开始构建集合可能更有效：

import random

def gen_2d(iterations, Chr_size, num_mut):
    randbps = set()
    while len(randbps) < iterations:
        listed = set()
        while len(listed) < num_mut:
            listed.add(random.choice(range(Chr_size)))
        randbps.add(tuple(sorted(listed)))
    return np.array(list(randbps))

这个函数从一个空集开始，在范围（Chr_size）中生成一个数字并将该数字添加到集合中。由于集合的属性，它不能再次添加相同的数字。它对 randbps 也做同样的事情，因此 randbps 的每个元素也是唯一的。

对于 np.random.choice 与 gen_2d 的仅一次迭代：

iterations=1000
Chr_size=1000000
num_mut=500

%timeit np.random.choice(range(Chr_size),num_mut,replace=False)
10 loops, best of 3: 141 ms per loop

%timeit gen_2d(1, Chr_size, num_mut)
1000 loops, best of 3: 647 µs per loop

Answer 2

基于this solution 中使用的技巧，这是一种在随机元素数组上使用argsort/argpartition 来模拟numpy.random.choice without replacement 以将randbps 作为二维数组提供给我们的方法-

np.random.rand(iterations,Chr_size).argpartition(num_mut)[:,:num_mut]

运行时测试-

In [2]: def original_app(iterations,Chr_size,num_mut):
   ...:     randbps=[]
   ...:     for k in range(iterations):
   ...:         listed=np.random.choice(range(Chr_size),num_mut,replace=False)
   ...:         randbps.append(listed)
   ...:     return randbps
   ...: 

In [3]: # Input params (scaled down version of params listed in question)
   ...: iterations=100
   ...: Chr_size=100000
   ...: num=50
   ...: 

In [4]: %timeit original_app(iterations,Chr_size,num)
1 loops, best of 3: 1.53 s per loop

In [5]: %timeit np.random.rand(iterations,Chr_size).argpartition(num)[:,:num]
1 loops, best of 3: 424 ms per loop