[Bioperl-l] Bio::Tools::SeqStats count_codons
Govind Chandra
govind.chandra at bbsrc.ac.uk
Tue Sep 27 11:15:51 EDT 2005
The count_codons method in Bio::Tools::SeqStats seems to be case
sensitive. See script below. Is this intended? Or am I doing something
silly.
### script begin ###
use Bio::SeqIO;
use Bio::Seq;
use Bio::Tools::SeqStats;
$seqio=Bio::SeqIO->new('-file' => 'Ssc.fas',
'-format' => 'fasta');
$seqobj=$seqio->next_seq();
open(ORFS, "<ordered_all_orfs.out");
while(<ORFS>) {
chomp();
($strand, $start, $end, $nt_len, $aa_len)=split(/\s+/, $_);
#print("$strand $start $end $nt_len $aa_len\n");
if($strand==1) {
$cdsobj=$seqobj->trunc($start, $end);
}
elsif($strand==-1) {
$cdsobj=$seqobj->trunc($start, $end)->revcom();
}
$tseq=uc($cdsobj->seq()); # works
#$tseq=$cdsobj->seq(); # does not work. given and ambiguous count only
$codobj=Bio::Seq->new('-seq' => $tseq,
'-alphabet' => 'dna'
);
print(">$strand $aa_len\n",$cdsobj->seq(),"\n");
$ss=Bio::Tools::SeqStats->new(-seq => $codobj);
$refcount = $ss->count_codons();
foreach $codon (keys(%$refcount)) {
print("$codon $refcount->{$codon}\n");
}
}
### script end ###
Cheers
Govind
John Innes Centre
Norwich UK.
More information about the Bioperl-l
mailing list