6.5.1 ClustalW¶
ClustalW是一个非常流行的进行多序列比对的命令行程序(其还有一个图形化的版本称之为ClustalX)。Biopython的 Bio.Align.Applications 模块包含这一多序列比对程序的打包程序。
我们建议你在Python中使用ClustalW之前在命令行界面下手动使用ClustalW,这样能使你更清楚这一程序的参数。你会发现Biopython打包程序非常严格地遵循实际的命令行API:
>>>from Bio.Align.Applications import ClustalwCommandline
>>>help(ClustalwCommandline)
作为最简单的一个例子,你仅仅需要一个FASTA格式的序列文件作为输入,例如: opuntia.fasta (你可以在线或者在Biopython/Doc/examples文件夹中找到该序列)。 opuntia.fasta 包含着7个prickly-pear的DNA序列(来自仙人掌科)。
ClustalW在默认情况下会产生一个包括所有输入序列的序列比对以及一个由输入序列名字构成的指导树(guide tree)。例如,用上述文件作为输入,ClustalW将会输出 opuntia.aln 和 opuntia.dnd 两个文件:
>>>from Bio.Align.Applications import ClustalwCommandline
>>>cline = ClustalwCommandline("clustalw2", infile="opuntia.fasta")
>>>print(cline)
clustalw2 -infile=opuntia.fasta
注意这里我们给出的执行文件名是 clustalw2 ,这是ClustalW的第二个版本(第一个版本的文件名为 clustalw )。ClustalW的这两个版本具有相同的参数,并且在功能上也是一致的。
你可能会发现,尽管你安装了ClustalW,以上的命令行却无法正确运行。你可能会得到“command not found”的错误信息(尤其是在Windows上)。这往往是由于ClustalW的运行程序并不在系统的工作目录PATH下(一个包含着运行程序路径的环境变量)。你既可以修改PATH,使其包括ClustalW的运行程序(不同系统需要以不同的方式修改),或者你也可以直接指定程序的绝对路径。例如:
>>>import os
>>>from Bio.Align.Applications import ClustalwCommandline
>>>clustalw_exe = r"C:\Program Files\new clustal\clustalw2.exe"
>>>clustalw_cline = ClustalwCommandline(clustalw_exe, infile="opuntia.fasta")
>>>assert os.path.isfile(clustalw_exe), "Clustal W executable missing"
>>>stdout, stderr = clustalw_cline()
注意,Python中 \n 和 \t 会被解析为一个新行和制表空白(tab)。然而,如果你将一个小写的“r”放在字符串的前面,这一字符串就将保留原始状态,而不被解析。这种方式对于指定Windows风格的文件名来说是一种良好的习惯。
Biopython在内部使用较新的 subprocess 模块来实现打包程序,而不是 os.system() 和 os.popen* 。
现在,我们有必要去了解命令行工具是如何工作的。当你使用一个命令行时,它往往会在屏幕上输出一些内容。这一输出可以被保存或重定向。在系统输出中,有两种管道(pipe)来区分不同的输出信息–标准输出(standard output)包含正常的输出内容,标准错误(standard error)显示错误和调试信息。同时,系统也接受标准输入(standard input)。这也是命令行工具如何读取数据文件的。当程序运行结束以后,它往往会返回一个整数。一般返回值为0意味着程序正常结束。
当你使用Biopython打包程序来调用命令行工具的时候,它将会等待程序结束,并检查程序的返回值。如果返回值不为0,Biopython将会提示一个错误信息。Biopython打包程序将会输出两个字符串,标准输出和标准错误。
在ClustalW的例子中,当你使用程序时,所有重要的输出都被保存到输出文件中。所有打印在屏幕上的内容(通过 stdout or stderr)可以被忽略掉(假设它已经成功运行)。
当运行ClustalW的时候,我们所关心的往往是输出的序列比对文件和指导树文件。ClustalW会自动根据输入数据的文件名来命名输出文件。在本例中,输出文件将是 opuntia.aln 。当你成功运行完ClustalW以后,你可以使用 Bio.AlignIO 来读取输出结果:
>>>from Bio import AlignIO
>>>align = AlignIO.read("opuntia.aln", "clustal")
>>>print(align)
SingleLetterAlphabet() alignment with 7 rows and 906 columns
TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273285|gb|AF191659.1|AF191
TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273284|gb|AF191658.1|AF191
TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273287|gb|AF191661.1|AF191
TATACATAAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273286|gb|AF191660.1|AF191
TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273290|gb|AF191664.1|AF191
TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273289|gb|AF191663.1|AF191
TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA gi|6273291|gb|AF191665.1|AF191
另一个输出文件 opuntia.dnd 中包含有一个newick格式的指导树,你可以使用Biopython中的 Bio.Phylo 来读取它:
>>>from Bio import Phylo
>>>tree = Phylo.read("opuntia.dnd", "newick")
>>>Phylo.draw_ascii(tree)
_______________ gi|6273291|gb|AF191665.1|AF191665
__________________________|
| | ______ gi|6273290|gb|AF191664.1|AF191664
| |__|
| |_____ gi|6273289|gb|AF191663.1|AF191663
|
_|_________________ gi|6273287|gb|AF191661.1|AF191661
|
|__________ gi|6273286|gb|AF191660.1|AF191660
|
| __ gi|6273285|gb|AF191659.1|AF191659
|___|
| gi|6273284|gb|AF191658.1|AF191658
13 章中详细介绍了如何使用Biopython对进化树数据进行处理。