view release on metacpan or  search on metacpan

lib/Bio/CUA/CUB/Calculator.pm  view on Meta::CPAN

			# Unknown parameter '$tag', ignored
		}
	}

	$self->no_atg(1); # exclude ATG in tAI calculation
	# check the input values
	# 1. make sure all the sense codons have CAI or tAI values if
	# provided

	return $self;
}

=head1 sequence input

all the following methods accept one of the following formats as
sequence input

=over

=item 1
 
 string of nucleotide sequence with length of 3N, 

=item 2

 sequence object which has a method I<seq> to get the sequence string,

=item 3

   a sequence file in fasta format

=item 4

   reference to a codon count hash, like
   $codons = { 
	   AGC => 50, 
       GTC => 124,
	   ...    ...
	   }.

=back

=head2 cai

 Title   : cai
 Usage   : $caiValue = $self->cai($seq);
 Function: calculate the CAI value for the sequence
 Returns : a number, or undef if failed
 Args    : see L</"sequence input">
 Note: codons without synonymous competitors are excluded in
 calculation.

=cut

sub cai
{
	my ($self, $seq) = @_;
	$self->_xai($seq, 'CAI');
}

# the real calculator of tAI or CAI as both have the same formula
sub _xai
{
	my ($self, $seq, $type) = @_;

	my $name;
	my $xaiHash;
	if($type =~ /cai/i)
	{
		$name = 'CAI';
		$xaiHash = $self->{"_cai_values"};
	}elsif($type =~ /tai/i)
	{
		$name = 'tAI';
		$xaiHash = $self->{"_tai_values"};
	}else
	{
		$self->throw("Unknown adaptation index type '$type'");
	}
	unless($xaiHash)
	{
		$self->warn("$name values for codons were not provided for",
		"this analyzer, so can not calculate $name for sequences");
		return undef;
	}

	my $codonList = $self->get_codon_list($seq) or return;

	my $xai = 0;
	my $seqLen = 0; # this excludes some unsuitable codons
	# get non-degenerative codons which are excluded in CAI
	my %nonDegCodons = map { $_ => 1 } $self->codons_by_degeneracy(1);
	my @senseCodons = $self->codon_table->all_sense_codons;
	foreach my $codon (@senseCodons)
	{
		next unless($codonList->{$codon}); # no observation of this codon
		# excluding non-degenerate codons for CAI calculation
		next if($nonDegCodons{$codon} and $type =~ /cai/i);
		unless(exists $xaiHash->{$codon})
		{
			$self->warn("Codon '$codon' is ignored")
			if($self->debug and ($self->no_atg? ($codon ne 'ATG') : 1));
			next;
		}
		my $cnt = $codonList->{$codon};
		# to overcome real number overflow, use log
		$xai += $cnt*log($xaiHash->{$codon});
		$seqLen += $cnt;
	}

	return undef unless($xai); # no codons with CAI/tAI

	$xai = exp($xai/$seqLen);
	return $xai;
}

=head2 fop

 Title   : fop
 Usage   : $fopValue = $self->fop($seq[,$withNonDegenerate]);
 Function: calculate the fraction of optimal codons in the sequence