[Biopython] Blast Two Sequences

[Peter Cock]

The NCBI BLAST+ arguments query and subject should be filenames, not
raw sequences.

The example you linked to on BioStars had an example by Ryan showing
how to first create two temp files which he called "seq1.fasta" and
"seq2.fasta" using Biopython's SeqIO.

Peter


On Thu, Mar 9, 2017 at 7:42 AM, Islam Amin <eng.islamamin at gmail.com> wrote:
> Dear All,
> I got a script from the following link : https://www.biostars.org/p/42687/
> about how to blast two sequences in biopython, I wonder if is it possible to
> pass the sequences direct to as following:
> NcbiblastpCommandline(query="FQTWEEFSRAAEKLYLADPMKVRV",subject="FQTWEEFSRAEKLYLADPMK",
> outfmt=5)()[0]
> unfortunately, I got the error when I'm trying to pass the sequence direct
> instead of read fasta file, is there any way to pass sequences direct.
>
>  the script is:
>
> 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
>
>
> _______________________________________________
> Biopython mailing list  -  Biopython at mailman.open-bio.org
> http://mailman.open-bio.org/mailman/listinfo/biopython