linux - biopython 肌肉命令行

Question

我正在使用 BioPython MuscleCommanLine 来比对子进程中的序列。肌肉的输入和输出是标准输入和标准输出。这行得通，但是一旦 popen 调用 muscle，我就会在屏幕上从 muscle 获得程序摘要。这大大减慢了程序速度，因为有数百万次调用子进程。

mcline = MuscleCommandline()
read_list = (SeqRecord(Seq(seq, IUPAC.unambiguous_dna), str(index)) for index, seq in enumerate(grouped_reads_list))

muscle = Popen(str(mcline), stdin=PIPE, stdout=PIPE, universal_newlines=True)

SeqIO.write(read_list, muscle.stdin, "fasta")  # Send sequences to Muscle in FASTA format.
muscle.stdin.close()

align = AlignIO.read(muscle.stdout, "fasta")  # Capture output from muscle and get it into FASTA format in an object.
print(align)
muscle.stdout.close()
exit("Testin Testing")

consensus_read = AlignInfo.SummaryInfo(align).dumb_consensus(threshold=0.6, ambiguous="N", consensus_alpha=IUPAC.ambiguous_dna)  # Create consensus from alignment object.

屏幕输出是

Robert C. Edgar 的 MUSCLE v3.8.31

http://www.drive5.com/muscle该软件捐赠给公共(public)领域。请引用:Edgar, R.C.核酸研究 32(5), 1792-97.

• 2 个序列，最大长度 133，平均长度 13300:00:00 10 MB(-1%) Iter 1 100.00% K-mer dist pass 100:00:00 10 MB(-1%) Iter 1 100.00% K-mer dist pass 200:00:00 12 MB(-1%) Iter 1 100.00% 对齐节点
  00:00:00 12 MB(-1%) Iter 1 100.00% 根对齐
• 6 个序列，最大长度 133，平均长度 133SingleLetterAlphabet() 对齐 2 行 133 列

Answer 1

我将其列为答案而不是编辑我的问题，因为有人可能会发现它有用。如果我犯了错误，请告诉我。问题似乎在于以这种方式使用 BioPython MuscleCommandLine 包装器。通过子进程调用时，我无法通过包装器将任何命令行选项传递给肌肉。我为此修改的代码如下。

cmd = ['muscle', "-quiet", "-maxiters", "1", "-diags"]

read_list = (SeqRecord(Seq(seq, IUPAC.unambiguous_dna), str(index)) for index, seq in enumerate(grouped_reads_list))

muscle = Popen(cmd, stdin=PIPE, stdout=PIPE, universal_newlines=True)

SeqIO.write(read_list, muscle.stdin, "fasta")  # Send sequences to Muscle in FASTA format.
muscle.stdin.close()

align = AlignIO.read(muscle.stdout, 'fasta')  # Capture output from muscle and get it into FASTA format in an object.

muscle.stdout.close()

consensus_read = AlignInfo.SummaryInfo(align).dumb_consensus(threshold=0.6, ambiguous="N", consensus_alpha=IUPAC.ambiguous_dna)
return str(consensus_read)