从FASTA文件提取序列到多个文件，文件基于单独fi中的头_id

Question

我正在寻找一个python解决方案，根据与单独文件中的头ID列表的匹配，将多个序列从FASTA文件提取到多个文件。在

这是发布在Extract sequences from a FASTA file based on entries in a separate file和{a2}上的问题的稍微复杂的版本，它们只为所有匹配项输出一个文件。

我是python新手，正在尝试找到一种方法：

• 获取一个包含将在fasta头文件中的字符串的文件
  
• 将所有与字符串匹配的记录写入单独的fasta文件
  

头文件如下：

CAP357_2030_09WPI、CAP357_2040_11WPI、CAP357_2050_13WPI等。。。在

我的fasta文件示例如下：

>；CAP357 U 2030 U 009wpi U v1v3_056_00002 U 000.4
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaaggg
>；CAP357_2040_011wpi_v1v3_008_00006_001.1
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt
>；CAP357_2040_011wpi_v1v3_030_00002_000.4
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt
>；CAP357_2040_011wpi_v1v3_004_00001_000.2
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt
>；CAP357_2050_013wpi_v1v3_047_00002_000.4
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt

预期输出

文件1:CAP357_2030_009wpi_v1v3.fasta
>；CAP357_2030_009wpi_v1v3_056_00002_000.4
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaaggg

文件2:CAP357_2040_011wpi_v1v3.fasta
>；CAP357_2040_011wpi_v1v3_008_00006_001.1
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt
>；CAP357_2040_011wpi_v1v3_030_00002_000.4
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt
>；CAP357_2040_011wpi_v1v3_004_00001_000.2
gtaaaattaaccacctgtgtcactctaaattgtacactgcaaagggt
等。。。在

这个代码来自上面的链接，但是我想要：
*写入单独外文件的匹配项
*如果可能，我不必单独指定每个输出文件（我最多有30个输出文件）

#!/usr/bin/env python

import sys
from Bio import SeqIO

input_file = sys.argv[1]
id_file = sys.argv[2]
output_file = sys.argv[3]

wanted = set(line.rstrip("\n").split(None,1)[0] for line in open(id_file))
print "Found %i unique identifiers in %s" % (len(wanted), id_file)

index = SeqIO.index(input_file, "fasta")
records = (index[r] for r in wanted)
count = SeqIO.write(records, output_file, "fasta")
assert count == len(wanted)

print "Saved %i records from %s to %s" % (count, input_file, output_file)

到目前为止，我已经想出了（下面的脚本），但是不知道如何手动指定所有的outfiles和变量（我在这里只包含了三个）

^{pr2}$

Answer 1

一些简短的建议：

如果所有标题都遵循相同的模式，则可以提取唯一的元素：

record.description.split("_")[1] 

（从“CAP357_2040_011wpi_v1v3_1_008_00006_001.1”得到“2040”）

如果使用dict，则可以收集记录集合：

^{pr2}$

然后，您可以将每个集合写入一个新文件：

file_name = "outfile%s" 
for (descr, records) in collected.items():   # iteritems in python2
    with open(os.path.join(file_path, file_name % descr), 'w') as f:
        SeqIO.write(records, f, 'fasta')

Answer 2

为了完整起见，以下是“最终”脚本：

#!/usr/bin/env python
# a script to extract fasta records from a fasta file to multiple separate fasta files based on a particular ID (time point) in a particular field, for a given delimiter
# to run, navigate to file location with command prompt and enter: python split_fasta_by_collections.py infile.fasta
from Bio import SeqIO
import os
import sys

records = SeqIO.parse(sys.argv[1], "fasta")
collected = {}
for record in records:
    descr = record.description.split("_")[1] # "_" sets the delimeter, "1" sets the field where counting starts at 0 for the first field
    try:
    collected[descr].append(record)
    except KeyError:
    collected[descr] = [record ,]

file_name = "outfile%s.fasta" 
file_path = os.getcwd() #sets the output file path to your current working directory

for (descr, records) in collected.items():  
    with open(os.path.join(file_path, file_name % descr), 'w') as f:
    SeqIO.write(records, f, 'fasta')