python - 使用 Python 进行高效重采样

Question

我想按如下方式重新采样序列:

fastadict = {"seq1" : "ATGCAGTCACGT", "seq2" : "ATGTGTGTACG"}

我写了下面的函数:

import sys
import random

def resampling_f(fastadict, seq, num):
    fastadict[seq] = fastadict[seq].replace("N","").replace("n","")
    l = []
    new_seq = ''.join([random.choice(fastadict[seq]) for i in range(num)]) 
    l.append(new_seq)
    return l

# Run function for 20 replicates:
for i in range(20):
    print resampling_f(fastadict, "seq1", 10)

这适用于示例中的小序列。在我的工作中，我需要对大约 100 万个字母(DNA 碱基，ACTG)进行 10000 次采样。这个功能对于这个目的来说太慢了。有没有一种更快的方法可以通过替换为 python 来获取采样？

Answer 1

使用 numpy.random 模块，它提供矢量化采样和 View 转换:

>>> n, k = 20, 10
>>> 
>>> np.random.choice((*fastadict['seq1'],), replace=True, size=(n * k,)).view(f'U{k}')
array(['GCAGAATGCT', 'GGAGGTGCAT', 'CACCATCATT', 'CGTGGTGTAC',
       'AGAATATCGG', 'GATTTTGGCC', 'GAACATAAGC', 'CGGGCCAAGC',
       'GTTGGTGTTT', 'GACCAATAAC', 'ACCCGTAGCC', 'GAATTCCCGG',
       'AACAGGTCAC', 'AGACAAGCAC', 'CACACTTGCC', 'CGTTTGTAAT',
       'CTAGCCCTCG', 'CTCGACATAT', 'GATGATTAGA', 'TCTATCCTCA'],
      dtype='<U10')

Python 2 版本:

>>> np.random.choice(tuple(fastadict['seq1']), replace=True, size=(n * k,)).view('S{k}'.format(k=k))

速度:

>>> from time import perf_counter
>>> n, k = 100, 1000000
>>> t0 = perf_counter(); x = np.random.choice((*fastadict['seq1'],), replace=True, size=(n * k,)).view(f'U{k}');t1 = perf_counter()
>>> t1-t0
1.29188625497045