选择汉明距离为零的读数

Question

我有一个 fastq 文件，说 reads.fastq。我有一个 7-mer 字符串列表。对于 reads.fastq 中的每个读取，我想检查它是否至少包含列表中的 7-mer 字符串之一。条件是，如果找到匹配项 (hamming distance ==0)，则将读取内容写入数组 chosen_reads，并匹配来自 fastq 文件的下一个读取内容。如果未找到匹配项，则循环继续，直到找到匹配项为止。输出数组由唯一读取组成，因为一旦找到第一个匹配项，匹配循环就会终止。我编写了以下代码，但输出数组中的读取不是唯一的，因为报告了所有汉明距离为零的匹配项。请提出修改建议：

def hamming(s1, s2):
    #Return the Hamming distance between equal-length sequences
    if len(s1) != len(s2):
        raise ValueError("Undefined for sequences of unequal length")

    return sum(ch1 != ch2 for ch1, ch2 in zip(s1, s2))

for x in Bio.SeqIO.parse("reads.fastq","fastq"):
        reads_array.append(x)

nmer = 7
l_chosen = ['gttattt','attattt','tgctagt']

chosen_reads = []
for x in reads_array:
    s2 = str(x.seq)
    for s in [s2[i:i+nmer] for i in range(len(s2)-nmer-1)]:
        for ds in l_chosen:    
            dist = hamming(ds,s)
            if dist == 0:
                print s2, s,ds,dist       
                chosen_reads.append(x)

Answer 1

当找到汉明距离为 0 的字符串时，您当前的代码不会跳出循环以从 reads.fastq 读取下一个 read，您应该使用标志来决定何时跳出out ，并在需要突破时为该标志分配 True 值 -

def hamming(s1, s2):
    #Return the Hamming distance between equal-length sequences
    if len(s1) != len(s2):
        raise ValueError("Undefined for sequences of unequal length")
    return sum(ch1 != ch2 for ch1, ch2 in zip(s1, s2))

for x in Bio.SeqIO.parse("reads.fastq","fastq"):
        reads_array.append(x)

nmer = 7

l_chosen = ['gttattt','attattt','tgctagt']
chosen_reads = []

for x in reads_array:
        s2 = str(x.seq)
        breakFlag = False
        for s in [s2[i:i+nmer] for i in range(len(s2)-nmer-1)]:
                for ds in l_chosen:
                        dist = hamming(ds,s)
                        if dist == 0:
                                print s2, s,ds,dist
                                chosen_reads.append(x)
                                breakFlag = True
                                break;
                if breakFlag:
                        break;

你确定要将 x 附加到 chosen_reads 中吗，这似乎是错误的，为了获得独特的匹配，也许你应该附加 s2 字符串和匹配 ds 而不是对吗？如果那是你想要的，你可以像下面那样将一个元组附加到 chosen_reads 而不是你当前的附加逻辑 -

chosen_reads.append((ds, s2))

Answer 2

如果我明白你在问什么，汉明距离就是试图准确地找到 3 个 "chosen" 字符串中的至少一个。像您一样进行迭代很慢，并且试图突破可能很丑陋。

我可能会建议 regex 在这里会有帮助。您可以自动创建匹配字符串：

import re
chosen_re = re.compile('|'.join(l_chosen))

chosen_reads = [x for x in reads_array if chosen_re.search(str(s.seq))]

你将很难超越正则表达式引擎的速度

选择汉明距离为零的读数

choosing reads with Hamming distance zero

python

bioinformatics

biopython

python-2.7