查找序列中不匹配的 DNA 条形码答案

【问题标题】：find a DNA barcode with mismatches in sequence查找序列中不匹配的 DNA 条形码
【发布时间】：2013-04-10 16:42:52
【问题描述】：

我有 36-nt 读取如下：atcttgttcaatggccgatcXXXXgtcgacaatcaa 在 fastq 文件中 XXXX 是不同的条形码。我想在文件中的确切位置（21 到 24）搜索条形码，并打印序列中最多 3 个不匹配的序列而不是条形码。

例如：我有条形码：aacg 在 fastq 文件中搜索位置 21 到 24 之间的条形码，允许序列中出现 3 个不匹配，例如：

atcttgttcaatggccgatcaacggtcgacaatcac # it has 1 mismatch
ttcttgttcaatggccgatcaacggtcgacaatcac # it has 2 mismatch
tccttgttcaatggccgatcaacggtcgacaatcac # it has 3 mismatch

我试图首先使用 awk 查找独特的行并查找不匹配的内容，但查找和查找它们对我来说非常乏味。

awk 'NR%4==2' 1.fq |sort|uniq -c|awk '{print $1"\t"$2}' > out1.txt

有什么快速的方法可以找到吗？

谢谢。

【问题讨论】：

我很困惑。条形码与核苷酸序列有什么关系？
最初我正在寻找特定位置的条形码，但我的计数非常低，并且序列中有 1 个不匹配我得到了高计数。所以，如果我在序列中给出不匹配，我会得到更多的序列（我想尝试最多 3 个）
所以您正在扫描barcodes？比如超市收银员用来识别商品价格的黑白条纹图案？因为我仍然不知道如何从条形码中获取 DNA。
@kevin en.wikipedia.org/wiki/DNA_barcoding

标签： python awk sed dna-sequence fastq

【解决方案1】：

使用 Python：

import re
seq="atcttgttcaatggccgatcaacggtcgacaatcaa"
D = [ c for c in seq ]
with open("input") as f:
    for line in f:
        line=line.rstrip('\n')
        if re.match(".{20}aacg", line):
            cnt = sum([ 1 for c,d in zip(line,D) if c != d])
            if cnt < 4:
                print cnt, line

【讨论】：

我只是看了一些输出 1 ATCTTGTTCAATGGCCGATCAACGGCCGACAATCAA 2 ATCTTGTTCAATGGCCGATCAACGGTCGACAATCAA 但第二个与序列相似，它说我有 2 个不匹配，为什么？
序列中必须有2个不匹配
original ATCTTGTTCAATGGCCGATCAACGGTCGACAATCAA 0 AGCTTGTTCAATGGCCGATCAACGGCCGACAATCAA 1 AGCTTGTTCAATGGCCGATCAACGGTCGACAATCAA 2 ATCTTGTTCAATGGCCGATCAACGGTCGACAATCAA 3 ATCTTGATCAATGGCCGATCAACGGTCGACAATCAA 这里是我得到@perreal的一些行
0, 2, 1 ,0, 1 是我得到的数据
是的。我得到了同样的东西。我只是想如果它找到相似的 0,1 错配序列与 1,2 错配序列与 2 和 3 错配序列与 3 @perreal

【解决方案2】：

使用 Python：

strs = "atcttgttcaatggccgatcaacggtcgacaatcaa"

with open("1.fq") as f:
    for line in f:
        if line[20:24] == "aacg":
            line = line.strip()
            mismatches = sum(x!=y for x, y  in zip(strs, line))
            if mismatches <= 3:
                print line, mismatches

atcttgttcaatggccgatcaacggtcgacaatcac 1
ttcttgttcaatggccgatcaacggtcgacaatcac 2
tccttgttcaatggccgatcaacggtcgacaatcac 3

【讨论】：

ATCTTGTTCAATGGCCGATCAACGGTCGACAATCAA ATCTTGTTCAATGGCCGATCAACGCTCGACAATCAA 3 AGCTTGTTCAATGGCCGATCAACGCTCGACAATCAA 2 我对输出感到困惑，这里有一些行，第一个是原始序列，数字 3 有 1 个不匹配，2 有 2 个不匹配。我想要的是它应该打印相似的序列（比如 0）和 1,2,3 不匹配的序列@Ashwini Chaudhary

【解决方案3】：

使用 python regex 模块可以指定不匹配的数量

import regex #intended as a replacement for re
from Bio import SeqIO
import collections

d = collections.defaultdict(list)
motif = r'((atcttgttcaatggccgatc)(....)(gtcgacaatcaa)){e<4}' #e<4 = less than 4 errors
records = list(SeqIO.parse(open(infile), "fastq"))
for record in records:
    seq = str(record.seq)
    match = regex.search(motif, seq, regex.BESTMATCH)
    barcode = match.group(3)
    sequence = match.group(0)
    d[barcode].append(sequence) # store as a dictionary key = barcode, value = list of sequences
for k, v in d.items():
    print("barcode = %s" % (k))
    for i in v:
        print("sequence = %s" % (i))

使用捕获组，第四组 (3) 将是条形码

【讨论】：