使用 Python 删除 FASTA 中的重复序列

Remove duplicated sequences in FASTA with Python

如果之前有人问过这个问题,我深表歉意,但我已经搜索了几天,在 Python 中找不到解决方案。

我有一个很大的 fasta 文件,包含 headers 和序列。

>cavPor3_rmsk_tRNA-Leu-TTA(m) range=chrM:2643-2717 5'pad=0 3'pad=0 strand=+ repeatMasking=none
GTTAAGGTGGCAGAGCCGGTAATTGCATAAAATTTAAGACTTTACTCTCA
GAGGTTCAACTCCTCTCCTTAACAC

>cavPor3_rmsk_tRNA-Gln-CAA_ range=chrM:3745-3815 5'pad=0 3'pad=0 strand=- repeatMasking=none
AGAGGGTCATAAAGGTTATGGGGTTGGCTTGAAACCAGCTTTAGGGGGTT
CAATTCCTTCCTCTCT

>cavPor3_rmsk_tRNA-Ser-TCA(m) range=chrM:6875-6940 5'pad=0 3'pad=0 strand=- repeatMasking=none
AGAGGGTCATAAAGGTTATGGGGTTGGCTTGAAACCAGCTTTAGGGGGTT
CAATTCCTTCCTCTCT

这是文件外观的一小部分。我只想保留第一个条目(header 和序列),如果您可以看到最后两个条目的序列相同。

输出将如下所示:

>cavPor3_rmsk_tRNA-Leu-TTA(m) range=chrM:2643-2717 5'pad=0 3'pad=0 strand=+ repeatMasking=none
GTTAAGGTGGCAGAGCCGGTAATTGCATAAAATTTAAGACTTTACTCTCA
GAGGTTCAACTCCTCTCCTTAACAC

>cavPor3_rmsk_tRNA-Gln-CAA_ range=chrM:3745-3815 5'pad=0 3'pad=0 strand=- repeatMasking=none
AGAGGGTCATAAAGGTTATGGGGTTGGCTTGAAACCAGCTTTAGGGGGTT
CAATTCCTTCCTCTCT

问题是 FASTA 文件的大小超过了 1 GB。我找到了通过基于重复 ID 删除重复项或使用 bash 来解决此问题的方法,但遗憾的是我无法在我的计算机上执行此操作。 此任务用于研究项目,而不是家庭作业或任务。

提前感谢您的帮助!

如果所有行都具有相同的格式,那么在序列之前有 6 个 space 分隔的字段,那么这很容易。您必须将所有唯一值存储在内存中。

memory = set()
for line in open('x.txt'):
    if len(line) < 5:
        continue
    parts = line.split(' ', 6)
    if parts[-1] not in memory:
        print( line.strip() )
        memory.add( parts[-1] )

这是从这里复制的:Remove Redundant Sequences from FASTA file in Python

使用 Biopython,但使用 fasta 文件,其中 headers 属于:

'> header' 类型见 FAsta Format Wiki

from Bio import SeqIO
import time

start = time.time() 

seen = []
records = []

for record in SeqIO.parse("INPUT-FILE", "fasta"):  
    if str(record.seq) not in seen:
        seen.append(str(record.seq))
        records.append(record)


#writing to a fasta file
SeqIO.write(records, "OUTPUT-FILE", "fasta")
end = time.time()

print(f"Run time is {(end- start)/60}") 


按照 MattMDo 的建议使用集合而不是列表更快:

seen = set()
records = []

for record in SeqIO.parse("b4r2.fasta", "fasta"):  
    if record.seq not in seen:
        seen.add(record.seq)
        records.append(record)

我有一个更长的使用 argparser 的,但速度较慢,因为计算序列可以 post 如果需要

如果您想在删除重复项的同时保留两者 header,您可以使用:

input1 = open("fasta.fasta")
dict_fasta = {record.id:record.seq for record in SeqIO.parse(input1,"fasta")}

tmp = {}
for key, value in dict_fasta.items():
  if value in tmp:
    tmp[value].append(key)
  else:
    tmp[value] = [ key ]