用python数一数fasta中的20毫秒数

Question

一个普通的fasta文件，读取长度为120 nt:'单个_地图.fa'

一个CSV文件包含10000个20 mer和每20个mer的计数：“20frequency”20梅斯.txt'像这样：

AAAAAGTATAGGAGATAGAA    35
AAAAATAGGAGGACTATTCA    26
AAAAATAGGAGGACTATTTA    24
AAAAATAGGAGGCCTATTCA    62

我想看完单曲_地图.fa，计算20frequency_20中所有20个MER的累计计数梅斯.txt对于每次读取，也就是说，对于read：

我想要61（35+26）

我的代码：

^{pr2}$

当我单独运行它们时，每个循环都能正常工作，但是没有像上面那样工作，有人能告诉我错误是从哪里来的吗？或者有没有更聪明、更有效的方法来做到这一点？提前谢谢！在

Answer 1

使用@MartinEvans dictionary的替代实现（不一定更好也更快），但使用re.findall()生成kmer进行测试，并使用map和{}代替（显式）内部循环：

from Bio import SeqIO
from re import findall
from itertools import repeat

kmers = {}

with open('20frequent_20mers.txt') as f_kmers:
    for line in f_kmers:
        kmer, count = line.strip().split('\t')
        kmers[kmer] = int(count)

for seq_record in SeqIO.parse("single_mapped.fa", "fasta"):
    print(seq_record.id)

    # use forward lookahead to make findall() find overlapping results;

    score_fre = sum(map(kmers.get, findall(r'(?=([ACTG]{20}))', str(seq_record.seq)), repeat(0)))

    print(score_fre)

Answer 2

使用现有的代码，您需要从一开始为每个序列和i值重新读取kmer文件。这将是非常缓慢的，应该避免。由于不将文件指针移回起始位置，因此它只能工作一次。在

可以通过在for row in kmer_list:行之前添加来移动文件指针：

file2.seek(0)

更好的方法是首先将所有kmer条目连同相应的计数一起加载到字典中。这样就可以快速查找它们：

^{pr2}$

如果在字典中找不到seq，则返回0的默认值。在