[Bioperl-l] Regarding Organism based search in Remote blast
Chris Fields
cjfields at illinois.edu
Fri Dec 4 09:59:17 EST 2009
Roopa,
At one point a couple of parameters differed between NCBI's web interface and our RemoteBlast-based BLAST interface to URLAPI (this should be indicated in your BLAST reports). See here:
http://thread.gmane.org/gmane.comp.lang.perl.bio.general/14155
Also, are the returned hits specific for the genome? You shoudl double-check; in some cases you have to set both HEADER and RETRIEVALHEADER to get the expected results (not sure why):
http://article.gmane.org/gmane.comp.lang.perl.bio.general/18737/match=remoteblast+ncbi
chris
On Dec 4, 2009, at 7:57 AM, Roopa Raghuveer wrote:
> Hello all,
>
> I am working on Remote blast.Here,I am trying to get 2 parameters into the
> remote blast code.They are
>
> 1.The input sequence that has to be sent to blast
>
> 2.Organism (The organism which has to be searched for ex:-Trypanasoma brucei
> etc.,)
>
> When I tried to take the organism parameter as an input from the
> user,through a web page,the Remote blast was not giving any results i.e., it
> says that there are no alignments found.
>
> But,when I hard coded the organism in the code,it gives me the results i.e.,
> 3hits.
>
> I could not understand this problem.Could any body please help me in this
> regard?
>
> My code is
>
> sub blastcode
> {
>
> $input1= $_[0];
>
> $organ= $_[1];
>
> open(NUC,'>',$nuc);
> print NUC $input1;
> close(NUC);
>
> my $prog = 'blastn';
> my $db = 'refseq_rna';
> my $e_val= '1e-10';
> my $organism= $organ;
>
> $gb = new Bio::DB::GenBank;
>
> my @params = ( '-prog' => $prog,
> '-data' => $db,
> '-expect' => $e_val,
> '-readmethod' => 'SearchIO',
> '-Organism' => $organism );
>
> open(OUTFILE,'>',$debugfile);
> print OUTFILE @params;
> close(OUTFILE);
>
>
> my $factory = Bio::Tools::Run::RemoteBlast->new(@params);
>
> #change a paramter
> $Bio::Tools::Run::RemoteBlast::HEADER{'ENTREZ_QUERY'} = '$organism[ORGN]';
> #change a paramter
> # $Bio::Tools::Run::RemoteBlast::HEADER{'ENTREZ_QUERY'} = '$input2[ORGN]';
>
> my $v = 1;
> #$v is just to turn on and off the messages
>
> my $str = Bio::SeqIO->new(-file => $nuc , '-format' => 'fasta' ,
> '-Organism' => $organism );
>
> while (my $input = $str->next_seq())
>
> {
> #Blast a sequence against a database:
> #Alternatively, you could pass in a file with many
> #sequences rather than loop through sequence one at a time
> #Remove the loop starting 'while (my $input = $str->next_seq())'
> #and swap the two lines below for an example of that.
>
> my $r = $factory->submit_blast($input);
>
> # my $r = $factory->submit_blast('amino.fa');
>
> print STDERR "waiting...." if($v>0);
>
> while ( my @rids = $factory->each_rid ) {
>
> foreach my $rid ( @rids ) {
>
> my $rc = $factory->retrieve_blast($rid);
>
> if( !ref($rc) )
> {
> if( $rc < 0 )
> {
> $factory->remove_rid($rid);
> }
> print STDERR "." if ( $v > 0 );
> sleep 5;
> }
> else {
> my $result = $rc->next_result();
> #save the output
> $blastdebugfile = $serverpath."/blastdebug_".time().".txt";
>
> # open(BLASTDEBUGFILE,'>',$debugfile);
> # print BLASTDEBUGFILE $result->next_hit();
> # close(BLASTDEBUGFILE);
>
> my $filename =
> $serverpath."/blastdata_".time().$result->query_name()."\.out";
>
> # open(DEBUGFILE,'>',$debugfile);
> # open(new,'>',$filename);
> # @arra=<new>;
> # print DEBUGFILE @arra;
> # close(DEBUGFILE);
> # close(new);
> $factory->save_output($filename);
>
> # open(BLASTDEBUGFILE,'>',$debugfile);
> # print BLASTDEBUGFILE "Hello $rid";
> # close(BLASTDEBUGFILE);
>
> $factory->remove_rid($rid);
>
> open(BLASTDEBUGFILE,'>',$blastdebugfile);
> print BLASTDEBUGFILE $organism;
> close(BLASTDEBUGFILE);
>
> # open(OUTFILE,'>',$outfile);
> # print OUTFILE "Test2 $result->database_name()";
> # close(OUTFILE);
>
> #$hit = $result->next_hit;
> #open(new,'>',$debugfile);
> #print $hit;
> #close(new);
>
> while ( my $hit = $result->next_hit ) {
>
> next unless ( $v > 0);
>
> # open(OUTFILE,'>',$debugfile);
> # print OUTFILE "$hit in while hits";
> # close(OUTFILE);
>
> my $sequ = $gb->get_Seq_by_version($hit->name);
> my $dna = $sequ->seq(); # get the sequence as a string
> push(@seqs,$dna);
> }
> }
> }
> }
> }
>
> #open(OUTFILE,'>',$debugfile);
> #print OUTFILE $seqs[0];
> #close(OUTFILE);
>
> return(@seqs);
> }
>
> Regards,
> Roopa.
> _______________________________________________
> Bioperl-l mailing list
> Bioperl-l at lists.open-bio.org
> http://lists.open-bio.org/mailman/listinfo/bioperl-l
More information about the Bioperl-l
mailing list