我有一个这样的fasta文件:test_fasta.fasta
>XXKHH_1
AAAAATTTCTGGGCCCC
>YYYXXKHH_1
TTAAAAATTTCTGGGCCCCGGGAAAAAA
>TTDTT_11
TTTGGGAATTAAACCCT
>ID_2SS
TTTGGGAATTAAACCCT
>YKHH_1
TTAAAAATTTCTGGGCCCCGGGAAAAAA
>YKHSH_1S
TTAAAAATTTCTGGGCCCCGGGAAAAAA
我想获得重复序列的计数并附加文件中每个序列的总计数(从最大到最小排序)并得到如下所示的结果:
>YYYXXKHH_1_counts3
TTAAAAATTTCTGGGCCCCGGGAAAAAA
>TTDTT_11_counts2
TTTGGGAATTAAACCCT
>XXKHH_1_counts1
AAAAATTTCTGGGCCCC
我有这个代码,它找到重复的序列并将它们连接在一起,但不是将它们连接在一起,我只想在id中附加重复项的计数,如上图所示。
from Bio import SeqIO
from collections import defaultdict
dedup_records = defaultdict(list)
for record in SeqIO.parse("test_fasta.fasta", "fasta"):
# Use the sequence as the key and then have a list of id's as the value
dedup_records[str(record.seq)].append(record.id)
with open("Output.fasta", 'w') as output:
for seq, ids in dedup_records.items():
# Join the ids and write them out as the fasta
output.write(">{}\n".format('|'.join(ids)))
output.write(seq + "\n")
由于您已经在输出循环的ids
列表中拥有每个重复记录的ID,因此您只需输出第一个ID(您显然想要的,根据您的预期输出),然后输出ids
列表的长度:
for seq, ids in sorted(dedup_records.items(), key=lambda t: len(t[1]), reverse=True):
output.write(">{}_counts{}\n".format(ids[0], len(ids)))
output.write(seq + "\n")