blast using python

Question

[santhoshhegde278]

hello, I am beginner of python I have two dna sequences in fasta format. I want to blast these two sequences.can anyone suggest me some solution in python. I have used NcbiblastnCommandline. but it is not worked.

from Bio.Blast.Applications import NcbiblastpCommandline
from StringIO import StringIO
from Bio.Blast import NCBIXML
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
from Bio import SeqIO

# Create two sequence files
seq1 = SeqRecord(Seq("FQTWEEFSRAAEKLYLADPMKVRVVLKYRHVDGNLCIKVTDDLVCLVYRTDQAQDVKKIEKF"),
                   id="seq1")
seq2 = SeqRecord(Seq("FQTWEEFSRAEKLYLADPMKVRVVLRYRHVDGNLCIKVTDDLICLVYRTDQAQDVKKIEKF"),
                   id="seq2")
SeqIO.write(seq1, "seq1.fasta", "fasta")
SeqIO.write(seq2, "seq2.fasta", "fasta")

# Run BLAST and parse the output as XML
output = NcbiblastpCommandline(query="seq1.fasta", subject="seq2.fasta", outfmt=5)()[0]
blast_result_record = NCBIXML.read(StringIO(output))

# Print some information on the result
for alignment in blast_result_record.alignments:
    for hsp in alignment.hsps:
        print '****Alignment****'
        print 'sequence:', alignment.title
        print 'length:', alignment.length
        print 'e value:', hsp.expect
        print hsp.query
        print hsp.match
        print hsp.sbjct

I found that error file not specified.

[st.ph.n]

If there is no particular reason to run BLAST through Biopython with XML output, I typically will run the blast command if I have downstream python code, using suprocess(), with tab-delimited output which is easier to parse IMO after with python.

[fishgolden]

Have you installed local blastp? Looking at your code, I assume you are using Windows.

Biopython tries to run just "blastp" not "blastp.exe" thus we must specify exact path to the blastp.exe (I think).

This is the line I ran blastp with your code.

output = NcbiblastpCommandline(cmd="C:\\apps\\ncbi-blast-2.6.0+\\bin\\blastp.exe",query="seq1.fasta", subject="seq2.fasta", outfmt=5)()[0]

C:\apps\ncbi-blast-2.6.0+\bin\blastp.exe is the executable of blastp I downloaded from NCBI FTP server.

[santhoshhegde278]

thank you very much.I will try this

[santhoshhegde278]

It is not working. I am getting same error.

[fishgolden]

Could you give me complete error message and your path of blastp.exe?

[User000]

yes, the code above works for me, so without the error message difficult to solve

Answer 1

[User000]

Why complicating it so much? if you have only two sequences just use NCBI blast web site, otherwise if you will have many sequences this is not the best method, you can run either as a command line or as a bash script.

#!/usr/bin/env bash
BLAST="/path/to/ncbi-blast-2.2.26+/bin"
REFERENCE="/path/to/reference.fa"
QUERY="/path/to/inputfile.fa"
OUTDIR="/path/to/output/directory"

echo "format..."
$BLAST/makeblastdb -in $REFERENCE -dbtype nucl
$BLAST/blastn -evalue 0.00001 -outfmt 6  -query $QUERY -db $REFERENCE -out $OUTDIR/outputfile

p.s. I also don't think parsing the .xml output of blast is a good idea

[Pierre Lindenbaum]

I also don't think parsing the .xml output of blast is a good idea

why ?

[User000]

on my personal experience I found it easier working with the tab delimited outputs. For a newbie in bioiformatics parsing an xml output might not be so straightforward, also if one has many sequences an xml format becomes heavier than tab. (See http://www.polarmicrobes.org/parsing-blast-xml-output/).