[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