1

我正在尝试编写一个代码,该代码将打开一个 fasta 文件并从不同的 fastq 文件中提取读取名称(标题)、序列(seq)和质量分数(qual),只有在 fasta 文件中找到它时,并写入将 fastq 信息放入一个新的 fastq 文件中。但是,我在如何编写最后一部分时遇到了麻烦(我在代码中遇到问题的地方加粗了)。可能有人知道如何编写这部分,或者我可以在哪里找到有关如何在 python 中输入的信息?

到目前为止,我有:

from sys import argv
from Bio.SeqIO.QualityIO import FastqGeneralIterator

script, merged_seqs, raw_seqs = argv
merged_from_raw = "merged_only.fastq"

merged_names = set() 
for line in open(merged_seqs): 
        if line[0] == ">":
                read_name = line.split()[0][1:] 
                merged_names.add(read_name) 

raw_fastq = raw_seqs

temp_handle = open(merged_from_raw, "w") 
for title, seq, qual in FastqGeneralIterator(open(raw_fastq)) :
        if title in merged_names:
                **handle.write() #this is where I don't know how to write what I need in python**
4

1 回答 1

1

除非您有特定的理由自己实现文件解析,否则最好使用 SeqIO 解析器来处理您的输入和输出文件。可能类似于以下内容(警告:我以前从未使用过 Bio,也没有测试过此代码):

from sys import argv
from Bio import SeqIO

output_filename = 'merged_only.fastq'
merged_seqs, raw_seqs = argv[1:2]

# Get fasta iterator, and read source fastq file into a dict-like object
merged_names = SeqIO.parse(merged_seqs, 'fasta')
source_seqs = SeqIO.index(raw_seqs, 'fastq')

filtered_seqs = (source_seqs[record.id] for record in merged_names if record.id in source_seqs)
SeqIO.write(filtered_seqs, output_filename, 'fastq')
于 2014-07-11T03:06:34.617 回答