使用 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 ]
如果之前有人问过这个问题,我深表歉意,但我已经搜索了几天,在 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 ]