[Bioperl-l] more genewise
Dimitar Kenanov
dimitark at bii.a-star.edu.sg
Fri May 7 05:03:58 UTC 2010
Hi guys,
another question about genewise. Is it possible to get the query seq and
the protein translation of the target seq somehow?
So, up to now i could not find a way to get the percent identity between
query and target(the protein translation) :( I spent some time on
CPAN and perldoc and even checked the code of several modules but still
no solution. Then i decided to extract the sequences out of the output
file and compare them somehow but i could not find a way and for that. I
found that the module 'Bio::Tools::Run::Genewise' is creating internal
temp output file which i cant access so i can parse it myself and
extract whatever.
Because with current implementation i cant access that temp output i
hacked a bit 'Bio::Tools::Run::Genewise' so i can pass my output file to
the constructor, like that:
my $factory = Bio::Tools::Run::Genewise->new( output => $tmpout);
#not "-output" cos the module currently doesnt like it
I modified the BEGIN section and the '_run' subroutine. My lines and the
originals are marked :
--------------
BEGIN {
@GENEWISE_PARAMS = qw( DYMEM CODON GENE CFREQ SPLICE GENESTATS INIT
SUBS INDEL INTRON NULL INSERT SPLICE_MAX_COLLAR
SPLICE_MIN_COLLAR
GW_EDGEQUERY GW_EDGETARGET GW_SPLICESPREAD
KBYTE HNAME ALG BLOCK DIVIDE GENER U V S T G E M);
@GENEWISE_SWITCHES = qw(HELP SILENT QUIET ERROROFFSTD TREV PSEUDO
NOSPLICE_GTAG
SPLICE_GTAG NOGWHSP GWHSP
TFOR TABS BOTH HMMER );
$OK_FIELD{OUTPUT}++; *#dimitar
* # Authorize attribute fields
foreach my $attr ( @GENEWISE_PARAMS, @GENEWISE_SWITCHES,
@OTHER_SWITCHES) { $OK_FIELD{$attr}++; }
}
-----------------------
-----------------------
my ($tfh1,$outfile1) = $self->io->tempfile(-dir=>$self->tempdir);
$self->debug("genewise command = $commandstring");
my $outfile2=$self->output; *#dimitar*
# my $status = system("$commandstring > $outfile1"); *#original*
my $status = system("$commandstring > $outfile2 "); *#dimitar*
$self->throw("Genewies call $commandstring crashed: $? \n") unless
$status==0;
# my $genewiseParser = Bio::Tools::Genewise->new(-file=> $outfile1);
*#original*
my $genewiseParser = Bio::Tools::Genewise->new(-file=> $outfile2);
*#dimitar*
-----------------------
More the method 'cds' from 'Bio::SeqFeature::Gene::Exon/I' gives nothing
back it doesnt matter what i tried. And i tried a lot :) Fortunately for
me i dont need that for now. But tried and didnt work so had to say.
Cheers
Dimitar
--
Dimitar Kenanov
Postdoctoral research fellow
Protein Sequence Analysis Group
Bioinformatics Institute
A*STAR, Singapore
tel: +65 6478 8514
More information about the Bioperl-l
mailing list