[Bioperl-l] Parsing Blast output for nohit sequences

Sammons, Scott zno6@cdc.gov
Wed, 31 Oct 2001 12:03:50 -0500


Greetings:

I am a relatively new bioperler.  I am trying to parse a large batch of
blast
searches and build a fasta file of the subjects that have no hits.  The code

below works fine when I want to do something with the hits.  But by the time

process_blast gets to the data, the nohits are already parsed out.  I can't 
seem to get to the blast parse while the nohits are still there. Any help
would 
be appreciated.

-Scott Sammons
Centers for Disease Control and Prevention
1600 Clifton Rd, NE
MS G36
Atlanta, GA 30333
404-639-3560 (Voice)
404-639-1331 (FAX)
ssammons@cdc.gov

==============
use Bio::Index::Fasta;
use Bio::SeqIO;
use Bio::Tools::Blast qw(:obj);

my ($inx, $seqobj, $outputfile, $blastobj, $params);

# parse blast output and create a fasta file from subject sequences with
nohits
# First, we create a index from the fasta file that was used to create the
blast output
# call this with  cat blastoutfile | nohit_blast2fasta.pl fastafile

sub process_blast {
        my $hit;
        my $blast_object = shift;
        if ($blast_object->num_hits('total') == 0) {
                $seqobj = $inx->fetch($blast_object->name);
                printf ">%s %s\n",$seqobj->display_id,$seqobj->desc;
                printf "%s\n", $seqobj->seq;
        }
        $blast_object->destroy;
}

$inx = Bio::Index::Fasta->new(
        -filename => "/tmp/t$$.idx",
        -write_flag => 1);
$inx->make_index(@ARGV);

$Blast->parse(  -parse => 1,
                -exec_func => \&process_blast,
                );