报告与 pairwise2 的最佳对齐

Question

我有一个读取的 fastq 文件，比方说 "reads.fastq"。我想将序列与保存为 fasta 文件 ref.faa 的字符串对齐。我为此使用以下代码

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)

但我只想报告每次读取的最佳比对。我如何从 a[2] 中的分数中选择这种对齐方式。我想选择具有最高 a[2]

值的对齐方式

Answer 1

您可以根据 a[2] 对路线进行排序：

for x in reads_array:
  alignments = pairwise2.align.globalms(
      str(refseq.seq).upper(), str(x.seq), 2, -1, -5, -0.05)
  sorted_alignments = sorted(alignments, key=operator.itemgetter(2))
        
  result.write(format_alignment(*sorted_alignments[0]))
  aligned_reads.append(x)

Answer 2

我知道这是一个老问题，但对于仍在寻找正确答案的任何人，请在对齐方法中添加 one_alignment_only=True 参数：

alignments =pairwise2.align.globalms(str(refseq.seq).upper(),
    str(x.seq),
    2,-1,-5,-0.05,
    one_alignment_only=True)

我必须在文档中进行一些挖掘才能找到它，但这会返回最佳分数。

报告与 pairwise2 的最佳对齐

reporting the best alignment with pairwise2

python

biopython