是否可以将字符串变量而不是文件传递给 BLAST 搜索?
Is it possible to pass a string variable to a BLAST search instead of a file?
我正在编写一个 python 脚本,如果可能的话,我想将查询序列信息作为字符串变量而不是 FASTA 格式文件传递到 blastn。
我使用 Biopython 的 SeqIO 将几个转录本名称存储为键,并将其序列存储为关联值。
所以它看起来像这样
transcripts = dict()
for record in SeqIO.parse("transcript_sequences.fasta", "fasta"):
transcripts[record.name] = record.seq
print transcripts
所以字典看起来像这样
{'var_F': Seq('CTTCATTCTCGTTTAGCGGCTGCTCGTGGAAATTTCGAAAAAATCTGAAACTAG...TGC', SingleLetterAlphabet())}
现在想把字典里面的序列信息解析成blast query和subject
subprocess.call("blastn -query " + transcript['var_F'] + ".txt" + " -subject " + transcript['var_B'] + " -outfmt 6 > tmp_blast.txt", shell=True)
我知道 blast 只接受文件名或文件位置字符串,但我只是想知道是否有解决方法。
提前感谢您阅读我的第一个问题:P
来自 BLAST 的 BioPython 模块:
http://biopython.org/DIST/docs/tutorial/Tutorial.html#htoc90
看来您是正确的,因为 biopython BLAST 不支持将 SeqIO 对象或生物序列作为 BLAST 函数调用的参数,或者当您使用 BLAST 二进制文件的 subprocess.call()
执行时。唯一接受的输入序列参数是文件名。来自教程:
>>> from Bio.Blast.Applications import NcbiblastxCommandline
>>> help(NcbiblastxCommandline)
...
>>> blastx_cline = NcbiblastxCommandline(query="opuntia.fasta", db="nr", evalue=0.001,
... outfmt=5, out="opuntia.xml")
>>> blastx_cline
NcbiblastxCommandline(cmd='blastx', out='opuntia.xml', outfmt=5, query='opuntia.fasta',
db='nr', evalue=0.001)
>>> print(blastx_cline)
blastx -out opuntia.xml -outfmt 5 -query opuntia.fasta -db nr -evalue 0.001
>>> stdout, stderr = blastx_cline()
因此,您唯一的选择是使用实际的 FASTA 文件作为输入。如果您想一次查询一个序列,则需要将每个序列保存到一个文件中。但是,我不建议这样做,除非您有理由这样做。我认为如果所有查询序列都在同一个文件中,BLAST 可能 执行得更快。您还可以使用 BioPython 读取生成的 BLAST 输出以迭代每个查询的结果,请参阅:
http://biopython.org/DIST/docs/tutorial/Tutorial.html#htoc92
上面的例子link:
如果您 运行 以其他方式进行 BLAST,并在文件 my_blast.xml 中具有 BLAST 输出(XML 格式),您需要做的就是打开阅读文件:
>>> result_handle = open("my_blast.xml")
>>> from Bio.Blast import NCBIXML
>>> blast_record = NCBIXML.read(result_handle)
或者,如果您有很多结果(即多个查询序列):
>>> from Bio.Blast import NCBIXML
>>> blast_records = NCBIXML.parse(result_handle)
就像Bio.SeqIO和Bio.AlignIO(见第5章和第6章)一样,我们有一对输入函数,read和parse,其中read是当你只有一个对象时,和parse 是一个迭代器,当你可以有很多对象时——但我们得到的不是 SeqRecord 或 MultipleSeqAlignment 对象,而是 BLAST 记录对象。
为了能够处理 BLAST 文件可能很大,包含数千个结果的情况,NCBIXML.parse() returns 一个迭代器。用简单的英语来说,迭代器允许您逐步执行 BLAST 输出,为每个 BLAST 搜索结果一条一条地检索 BLAST 记录:
>>> from Bio.Blast import NCBIXML
>>> blast_records = NCBIXML.parse(result_handle)
>>> blast_record = next(blast_records)
# ... do something with blast_record
>>> blast_record = next(blast_records)
# ... do something with blast_record
>>> blast_record = next(blast_records)
# ... do something with blast_record
>>> blast_record = next(blast_records)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
StopIteration
# No further records
您可以在 运行 blast 命令时将序列作为 fasta 格式的字符串传递到标准输入。
query_string = '>' + record.description + '\n' + str(record.seq)
blast_command = Bio.Blast.Applications.NcbiblastnCommandline(cmd='blastn', out=record.id + '.xml', outfmt=5, db=database_name, evalue=0.001)
stdout, stderr = blast_command(stdin=query_string)
我正在编写一个 python 脚本,如果可能的话,我想将查询序列信息作为字符串变量而不是 FASTA 格式文件传递到 blastn。
我使用 Biopython 的 SeqIO 将几个转录本名称存储为键,并将其序列存储为关联值。
所以它看起来像这样
transcripts = dict()
for record in SeqIO.parse("transcript_sequences.fasta", "fasta"):
transcripts[record.name] = record.seq
print transcripts
所以字典看起来像这样
{'var_F': Seq('CTTCATTCTCGTTTAGCGGCTGCTCGTGGAAATTTCGAAAAAATCTGAAACTAG...TGC', SingleLetterAlphabet())}
现在想把字典里面的序列信息解析成blast query和subject
subprocess.call("blastn -query " + transcript['var_F'] + ".txt" + " -subject " + transcript['var_B'] + " -outfmt 6 > tmp_blast.txt", shell=True)
我知道 blast 只接受文件名或文件位置字符串,但我只是想知道是否有解决方法。
提前感谢您阅读我的第一个问题:P
来自 BLAST 的 BioPython 模块:
http://biopython.org/DIST/docs/tutorial/Tutorial.html#htoc90
看来您是正确的,因为 biopython BLAST 不支持将 SeqIO 对象或生物序列作为 BLAST 函数调用的参数,或者当您使用 BLAST 二进制文件的 subprocess.call()
执行时。唯一接受的输入序列参数是文件名。来自教程:
>>> from Bio.Blast.Applications import NcbiblastxCommandline
>>> help(NcbiblastxCommandline)
...
>>> blastx_cline = NcbiblastxCommandline(query="opuntia.fasta", db="nr", evalue=0.001,
... outfmt=5, out="opuntia.xml")
>>> blastx_cline
NcbiblastxCommandline(cmd='blastx', out='opuntia.xml', outfmt=5, query='opuntia.fasta',
db='nr', evalue=0.001)
>>> print(blastx_cline)
blastx -out opuntia.xml -outfmt 5 -query opuntia.fasta -db nr -evalue 0.001
>>> stdout, stderr = blastx_cline()
因此,您唯一的选择是使用实际的 FASTA 文件作为输入。如果您想一次查询一个序列,则需要将每个序列保存到一个文件中。但是,我不建议这样做,除非您有理由这样做。我认为如果所有查询序列都在同一个文件中,BLAST 可能 执行得更快。您还可以使用 BioPython 读取生成的 BLAST 输出以迭代每个查询的结果,请参阅:
http://biopython.org/DIST/docs/tutorial/Tutorial.html#htoc92
上面的例子link:
如果您 运行 以其他方式进行 BLAST,并在文件 my_blast.xml 中具有 BLAST 输出(XML 格式),您需要做的就是打开阅读文件:
>>> result_handle = open("my_blast.xml")
>>> from Bio.Blast import NCBIXML
>>> blast_record = NCBIXML.read(result_handle)
或者,如果您有很多结果(即多个查询序列):
>>> from Bio.Blast import NCBIXML
>>> blast_records = NCBIXML.parse(result_handle)
就像Bio.SeqIO和Bio.AlignIO(见第5章和第6章)一样,我们有一对输入函数,read和parse,其中read是当你只有一个对象时,和parse 是一个迭代器,当你可以有很多对象时——但我们得到的不是 SeqRecord 或 MultipleSeqAlignment 对象,而是 BLAST 记录对象。
为了能够处理 BLAST 文件可能很大,包含数千个结果的情况,NCBIXML.parse() returns 一个迭代器。用简单的英语来说,迭代器允许您逐步执行 BLAST 输出,为每个 BLAST 搜索结果一条一条地检索 BLAST 记录:
>>> from Bio.Blast import NCBIXML
>>> blast_records = NCBIXML.parse(result_handle)
>>> blast_record = next(blast_records)
# ... do something with blast_record
>>> blast_record = next(blast_records)
# ... do something with blast_record
>>> blast_record = next(blast_records)
# ... do something with blast_record
>>> blast_record = next(blast_records)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
StopIteration
# No further records
您可以在 运行 blast 命令时将序列作为 fasta 格式的字符串传递到标准输入。
query_string = '>' + record.description + '\n' + str(record.seq)
blast_command = Bio.Blast.Applications.NcbiblastnCommandline(cmd='blastn', out=record.id + '.xml', outfmt=5, db=database_name, evalue=0.001)
stdout, stderr = blast_command(stdin=query_string)