[Biopython] Blast Two Sequences
Islam Amin
eng.islamamin at gmail.com
Thu Mar 9 23:15:30 UTC 2017
Many thanks Peter, I thought that It may be possible to pass raw sequences
instead of passing fasta file.
On Fri, Mar 10, 2017 at 9:43 AM, Peter Cock <p.j.a.cock at googlemail.com>
wrote:
> 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("FQTWEEFSRAAEKLYLADPMKVRVVLKYRH
> VDGNLCIKVTDDLVCLVYRTDQAQDVKKIEKF"),
> > id="seq1")
> > seq2 =
> > SeqRecord(Seq("FQTWEEFSRAEKLYLADPMKVRVVLRYRHV
> DGNLCIKVTDDLICLVYRTDQAQDVKKIEKF"),
> > 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
>
--
Best Regards,
Islam Amin.
www.egyptscience.net
Scientific Research Group in EGYPT
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.open-bio.org/pipermail/biopython/attachments/20170310/a0c01a42/attachment.html>
More information about the Biopython
mailing list