[Bioperl-l] Writing blast report with HitTableWriter

Jason Stajich jason@cgt.mc.duke.edu
Mon, 24 Jun 2002 12:04:45 -0400 (EDT)


Steve wrote the texttable dumper for his objects created by psiblast
parser (will parse both blast & psiblast) but the dumper expected methods
that were not in the Generic objects produced by SearchIO::blast.

However, this is *almost* fixed to work on general GenericXXX objects that
are created by Bio::SearchIO::blast - live code in CVS will have these
changes (main trunk only I think).  You can try everything using Steve's
parser - -format => 'psiblast' and it might work okay.



On Mon, 24 Jun 2002, Vilanova,David,LAUSANNE,NRC/BS wrote:

> Hello guys,
> I've been playing with a script i found on the exemple section. I try to
> customize it but i have a problem.
> Basically all the fields are equal to 0 in the output after the parsing of
> the blast report.
> Does anyone know why ???
> Thanks,
> David
>
> Output below:
> ---	-------	-------
> 	0		0	0	0.0e+00	0.00	0	0
>
> 	0		0	0	0.0e+00	0.00	0	0
>
> 	0		0	0	0.0e+00	0.00	0	0
>
>
> Script:
>
> use Bio::SeqIO;
> use Bio::Tools::Run::StandAloneBlast;
> use Bio::SearchIO;
> use Bio::SearchIO::Writer::HitTableWriter;
>
> #Check STDIN
> die "Usage: perl script.pl input_file_to_blast output_file\n" unless @ARGV
> eq '2';
>
> ($file,$output) = @ARGV;
> $raw_blast_output = "raw_blast.".$output;
>
> printf ("Preparing files....\n");
>
> open (BLAST_REPORT,">$output") or die "cannot open $output for writing\n";
>
>
> #Define blast parameters.
> my $blast = Bio::Tools::Run::StandAloneBlast->new(
> 						  database => 'nrdb',
> 						  program => 'blastx',
> 						  '_READMETHOD' => 'Blast',
> 						  # to automatically create
> a
> 						  # Bio::Tools::SearchIO
> object rather
> 						  # than a
> Bio::Tools::BPlite one
>
> 						 );
>
> # Here I add more blast parameters
> $blast->a(8); #nb of cpu
> $blast->e(0.0000000001); #e-10
> #$blast->m(7); #generate XML output
> $blast->o($raw_blast_output);
>
> #*********open file containing the sequences to blast**********
>
> my $Seq_in = Bio::SeqIO->new ( -file   => $file,
> 			       -format => 'fasta');
>
>
> $seq = $Seq_in->next_seq(); #Read sequence
> $blast_out = $blast->blastall($seq); #Run blast
>
>
>
> # These are the columns that will be in the output table of BLAST results.
> my @columns = qw(
> 		 query_name
> 		 query_length
>                  hit_name
>                  hit_length
> 		 num_hsps
>                  expect
>                  frac_identical_query
>                  length_aln_query
>                  gaps_total
>                  strand_query
>                  strand_hit
> 		);
>
>
>
>
>
> my $in     =  Bio::SearchIO->new( -format => 'blast',
> 				  -hit_filter => $filt_func,
> 				  -file => $raw_blast_output,
> 				 );
>
> my $writer = Bio::SearchIO::Writer::HitTableWriter->new( -columns =>
> \@columns,
> 						       );
> my $out    = Bio::SearchIO->new(  -writer => $writer,
> 				  -file   => ">$output" );
> my $hit_count = 0;
>
> while ( my $blast_result = $in->next_result() ) {
>   #printf STDERR "\nReport %d: $blast\n", $in->report_count;
>
>   if( $blast_result->hits ) {
>     $hit_count++;
>     $out->write_result($blast_result, $hit_count==1 );
>   }
>   else {
>     print STDERR "Hitless Blast Report: ";
>     print STDERR ($blast->no_hits_found ? "\n" : "(filtered)\n");
>   }
> }
>
> _______________________________________________
> Bioperl-l mailing list
> Bioperl-l@bioperl.org
> http://bioperl.org/mailman/listinfo/bioperl-l
>

-- 
Jason Stajich
Duke University
jason at cgt.mc.duke.edu