This is a test version of Biostars. For the public version, visit https://www.biostars.org.
reporting the best alignment(s) only in pairwise2

I have a fastq file of reads, say "reads.fastq". I want to align the sequences to a string saved as a fasta file ref.faa. I am using the following code for this

reads_array = []

for x in Bio.SeqIO.parse("reads.fastq","fastq"):
    reads_array.append(x)

for x in Bio.SeqIO.parse("ref.faa","fasta"):
    refseq = x

result = open("alignments_G10_final","w")

aligned_reads = []

for x in reads_array:
    alignments =pairwise2.align.globalms(str(refseq.seq).upper(),str(x.seq),2,-1,-5,-0.05)
    for a in alignments:
        result.write(format_alignment(*a))
        aligned_reads.append(x)

But I want to report only the best alignment for each read. How can I choose this alignment from the scores in a[2]. I want to report the best alignment(s) with the highest value of score

alignment next-gen-sequencing

1 answer

just add at the end: one_alignment_only=True

Log in to answer this question.