如何扩展模糊的dna序列

jrj*_*rjc 5 python dna-sequence biopython

假设您的DNA序列如下:

AATCRVTAA
Run Code Online (Sandbox Code Playgroud)

其中RV是DNA核苷酸的暧昧值,其中R任一代表AGV表示ACG

是否存在Biopython方法来生成上述歧义序列所代表的序列的所有不同组合?

例如,在这里,输出将是:

AATCAATAA
AATCACTAA
AATCAGTAA
AATCGATAA
AATCGCTAA
AATCGGTAA
Run Code Online (Sandbox Code Playgroud)

Jiv*_*van 5

也许更短更快的方法,因为无论如何这个函数将用于非常大的数据:

from Bio import Seq
from itertools import product

def extend_ambiguous_dna(seq):
   """return list of all possible sequences given an ambiguous DNA input"""
   d = Seq.IUPAC.IUPACData.ambiguous_dna_values
   return [ list(map("".join, product(*map(d.get, seq)))) ]
Run Code Online (Sandbox Code Playgroud)

Usingmap允许您在 C 中而不是在 Python 中执行循环。这应该比使用普通循环甚至列表推导式要快得多。

现场测试

用一个简单的 dict asd而不是返回的那个ambiguous_na_values

from itertools import product
import time

d = { "N": ["A", "G", "T", "C"], "R": ["C", "A", "T", "G"] }
seq = "RNRN"

# using list comprehensions
lst_start = time.time()
[ "".join(i) for i in product(*[ d[j] for j in seq ]) ]
lst_end = time.time()

# using map
map_start = time.time()
[ list(map("".join, product(*map(d.get, seq)))) ]
map_end = time.time()

lst_delay = (lst_end - lst_start) * 1000
map_delay = (map_end - map_start) * 1000

print("List delay: {} ms".format(round(lst_delay, 2)))
print("Map delay: {} ms".format(round(map_delay, 2)))
Run Code Online (Sandbox Code Playgroud)

输出:

# len(seq) = 2:
List delay: 0.02 ms
Map delay: 0.01 ms

# len(seq) = 3:
List delay: 0.04 ms
Map delay: 0.02 ms

# len(seq) = 4
List delay: 0.08 ms
Map delay: 0.06 ms

# len(seq) = 5
List delay: 0.43 ms
Map delay: 0.17 ms

# len(seq) = 10
List delay: 126.68 ms
Map delay: 77.15 ms

# len(seq) = 12
List delay: 1887.53 ms
Map delay: 1320.49 ms
Run Code Online (Sandbox Code Playgroud)

显然map更好,但只是提高了 2 或 3 倍。可以肯定的是,它可以进一步优化。