view release on metacpan or  search on metacpan

lib/Bio/Tools/CodonOptTable.pm  view on Meta::CPAN

    $alphabet = $seq_data->alphabet;

    return ( $seq, $id, $desc, $alphabet );
}

sub _read_localfile {
    my ( $self, $file, $format ) = @_;

    my ( $seq, $id, $alphabet );

    my $inputstream = Bio::SeqIO->new(
        -file   => $file,
        -format => $format
    );
    my $input = $inputstream->next_seq();

    $seq      = $input->seq();
    $id       = $input->id;
    $alphabet = $input->alphabet;

    return ( $seq, $id, $alphabet );
}

sub rscu_rac_table {
    my $self = shift;

    my $codons = $self->{'codons'};
    my ( @codons, $max_rscu ) = $self->_calculate_rscu();
    my $rscu_rac = $self->_calculate_rac( @codons, $max_rscu );
    my @sorted_codons_by_aa =
      sort { $a->{aa_name} cmp $b->{aa_name} } @$rscu_rac;

    return \@sorted_codons_by_aa;
}

sub _map_codon_iupac {
    my $self = shift;

    my $myCodonTable = Bio::Tools::CodonTable->new();
    my $codons       = $self->{'codons'};
    $self->{'genetic_code'} = 1 if ( !$self->{'genetic_code'} );

    $myCodonTable->id( $self->{'genetic_code'} );

    my $myCodons;

    foreach my $single_codon ( keys %$codons ) {
        my $aa_name_abri = $myCodonTable->translate($single_codon);
        my $aa_name      = $AMINOACID{$aa_name_abri};
        $myCodons->{$single_codon} = {
            'frequency' => $codons->{$single_codon},
            'aa_name'   => $aa_name,
        };
    }
    $self->{'codons'} = $myCodons;

    return 1;
}

#    Function : Calculate the RSCU(Relative Synonymous Codons Uses).
#    Note     : The formula is used in the following references.
#	 http://www.pubmedcentral.nih.gov/articlerender.fcgi?tool=pubmed&pubmedid=3547335

sub _calculate_rscu {
    my $self = shift;

    my $codons = $self->{'codons'};

    my ( $rscu, @myCodons, %rscu_max_table );

    foreach my $each_codon ( keys %$codons ) {
        my $amino = $codons->{$each_codon}->{'aa_name'};
        my $freq  = $codons->{$each_codon}->{'frequency'};
        my $count;
        my $all_freq_aa = 0;
        if ($amino) {
            foreach my $goforall ( keys %$codons ) {
                if ( $amino && ( $codons->{$goforall}->{'aa_name'} eq $amino ) )
                {
                    $all_freq_aa += $codons->{$goforall}->{'frequency'};
                    $count++;
                }
            }
            $rscu = $all_freq_aa > 0 ? $count * $freq / $all_freq_aa : 0.00;
        }
        if ( !defined( $rscu_max_table{$amino} )
            || ( $rscu_max_table{$amino} < $rscu ) )
        {
            $rscu_max_table{$amino} = $rscu;
        }
        push @myCodons,
          {
            codon         => $each_codon,
            frequency     => $freq,
            aa_name       => $amino,
            rscu          => $rscu,
            total_aa_comb => $count,
            all_fre_aa    => $all_freq_aa,
          };
    }
    return ( \@myCodons, \%rscu_max_table );
}

#    Function : Calculate the RAC (Relative Adaptiveness of a Codon).
#    Note     : The formula is used in the following references.
#	 http://www.pubmedcentral.nih.gov/articlerender.fcgi?tool=pubmed&pubmedid=3547335

sub _calculate_rac {
    my ( $self, $codons, $max_rscu ) = @_;
    my ( $rac, @myCodons );

    foreach my $each_codon (@$codons) {
        my $amino = $each_codon->{'aa_name'};
        my $rscu  = $each_codon->{'rscu'};
        if ($amino) {
            my $max = $max_rscu->{$amino};
            $rac = $max > 0 ? $rscu / $max : 0.00;
            push @myCodons,
              {
                'codon'     => $each_codon->{'codon'},
                'frequency' => $each_codon->{'frequency'},
                'aa_name'   => $amino,
                'rscu'      => sprintf( "%.2f", $rscu ),
                'rac'       => sprintf( "%.2f", $rac ),
              };
        }
    }
    return ( \@myCodons );
}

# The CAI index is defined as the geometric mean of these relative adaptiveness values.
sub calculate_cai {
    my ( $self, $rscu_rac ) = @_;

    my $count;
    my $w;
    foreach my $codon (@$rscu_rac) {
        $w += $codon->{rac};
        $count++;
    }
    my $cai = $w * ( 1 / $count );

# A CAI of 1.0 is considered to be perfect in the desired expression organism, and a
# CAI of >0.9 is regarded as very good, in terms of high gene expression level.
    return sprintf( "%.2f", $cai );
}

sub prefered_codon {
    my ( $self, $codons ) = @_;
    my $prefered_codon;
    for my $each_aa (@$codons) {
        my $aa_name   = $each_aa->{'aa_name'};
        my $rscu      = $each_aa->{'rscu'};
        my $codon     = $each_aa->{'codon'};
        my $frequency = $each_aa->{'frequency'};

        if ( !defined( $prefered_codon->{$aa_name} )
            || ( $prefered_codon->{$aa_name} < $rscu ) )
        {
            $prefered_codon->{$aa_name} = $codon;
        }
    }
    return $prefered_codon;
}