我正在解决一个DNA序列的模式匹配问题,汉明距离达d。Regex救了我一命。但现在我遇到了另一个问题。给定一个长的DNA序列,我必须找到最常见的错配k-mer,最多有d个错配。这里,k-mer指长度k的子序列
注意:DNA序列只能用四个字母来表示:{A,C,G,T}
比如说,
DNA sequence= "AGGC"
k = 3
d = 1
这里,只有两个k-mers是可能的:“AGG”、“GGC”
现在,我可以通过为“AGG”和“GGC”运行以下代码,将它们分别排列为1个不匹配项
def permute_one_nucleotide(motif, alphabet={"A", "C", "G", "T"}):
import itertools
return list(
set(
itertools.chain.from_iterable(
[
[
motif[:pos] + nucleotide + motif[pos + 1 :]
for nucleotide in alphabet
]
for pos in range(len(motif))
]
)
)
)
“AGG”将给出:
['TGG', 'ATG', 'AGG', 'GGG', 'AGT', 'CGG', 'AGC', 'AGA', 'ACG', 'AAG']
“GCC”将给出:
['GCC', 'GAC', 'GGT', 'GGA', 'AGC', 'GTC', 'TGC', 'CGC', 'GGG', 'GGC']
然后我可以使用Counter
查找最频繁的k-mer。但是,这只有在d = 1
时才有效。如何将其推广到任何d <= k
这是一种计算成本很高的方法。但是,是的,应该会得到你想要的。 我在这里做的是计算所有与hamming dist 1的不匹配。然后用ham dist 1从上一次失配计算新的失配,并递归到d
发现性能更好的代码,速度快了近10-20倍
CPU时间:用户1分钟32秒,系统1.81秒,总计1分钟34秒 墙壁时间:1分钟34秒
相关问题 更多 >
编程相关推荐