Skip to main content
 首页 » 编程设计

python之在 python 循环中从交替文件打印行

2025年05月04日22jpfss

我正在尝试使用 python 在两个单独的文件中找到感兴趣的四行 block ,然后按受控顺序打印出其中一些行。下面是两个输入文件和所需输出文件的示例。请注意,Input.fasta 中的 DNA 序列与 Input.fastq 中的 DNA 序列不同,因为 .fasta 文件已被读取更正。

输入.fasta

>read1 
AAAGGCTGT 
>read2 
AGTCTTTAT 
>read3 
CGTGCCGCT 

输入.fastq

@read1 
AAATGCTGT 
+ 
'(''%$')) 
@read2 
AGTCTCTAT 
+ 
&---+2010 
@read3 
AGTGTCGCT 
+ 
0-23;:677 

DesiredOutput.fastq

@read1 
AAAGGCTGT 
+ 
'(''%$')) 
@read2 
AGTCTTTAT 
+ 
&---+2010 
@read3 
CGTGCCGCT 
+ 
0-23;:677 

基本上我需要序列行“AAAGGCTGT”, 来自“input.fasta”的“AGTCTTTAT”和“CGTGCCGCT”以及来自“input.fastq”的所有其他行。这允许将质量信息恢复到读取更正的 .fasta 文件。

这是我最接近的失败尝试:

fastq = open(Input.fastq, "r") 
fasta = open(Input.fasta, "r") 
 
ReadIDs = [] 
IDs = [] 
 
with fastq as fq: 
    for line in fq: 
        if "read" in line:   
            ReadIDs.append(line) 
            print(line.strip()) 
            for ID in ReadIDs: 
                IDs.append(ID[1:6]) 
            with fasta as fa: 
                for line in fa: 
                    if any(string in line for string in IDs): 
                        print(next(fa).strip()) 
            next(fq) 
            print(next(fq).strip()) 
            print(next(fq).strip()) 

我想我遇到了麻烦,因为我试图在同一个循环中嵌套对两个不同文件的“with”调用。这会正确打印 read1 所需的行,但不会继续遍历其余行并抛出错误“ValueError:已关闭文件的 I/O 操作”

请您参考如下方法:

我建议你使用 Biopython ,这将为您省去很多麻烦,因为它为这些文件格式提供了很好的解析器,它不仅可以处理标准情况,还可以处理多行 fasta。

这是一个用相应的fasta序列行替换fastq序列行的实现:

from Bio import SeqIO 
 
fasta_dict = {record.id: record.seq for record in 
              SeqIO.parse('Input.fasta', 'fasta')} 
 
def yield_records(): 
    for record in SeqIO.parse('Input.fastq', 'fastq'): 
        record.seq = fasta_dict[record.id] 
        yield record 
 
SeqIO.write(yield_records(), 'DesiredOutput.fastq', 'fastq') 

如果你不想使用标题而只是依赖顺序那么解决方案更简单并且内存效率更高(只需确保顺序和记录数相同),无需定义字典首先,一起遍历记录:

fasta_records = SeqIO.parse('Input.fasta', 'fasta') 
fastq_records = SeqIO.parse('Input.fastq', 'fastq') 
 
def yield_records(): 
    for fasta_record, fastq_record in zip(fasta_records, fastq_records): 
        fastq_record.seq = fasta_record.seq 
        yield fastq_record