[Bioperl-l] Getting pairwise alignment scores for existing multiple alignment

Alexey Morozov alexeymorozov1991 at gmail.com
Wed Oct 1 01:12:59 UTC 2014


Thanks Chris,
the module you mentioned can only calculate score for nucleic acids, and I
am in fact working with aminoacids alignment. However, I've written the sub
for my needs (it's below, in case someone else needs it).

sub score_aln
        {
        #Takes two alignment lines AS STRINGS;
        #Assumes $matrix (Bio::Matrix::Generic), $open_penalty and
$extend_penalty to be declared as global variables

        die "Different lengths of sequences:\n$_[0]\n$_[1]\n" if
length$_[0]!=length$_[1];
        my @line1=split(//,$_[0]);
        my @line2=split(//,$_[0]);
        my $score=0;
        my $extend=0;
        #Below we count score for each pair (two AA, new gap or gap
extension)
        CHAR:for (my $j=0;$j<scalar(@line1);$j++)
                {
                if (($line1[$j] eq '-')or($line2[$j] eq '-'))
                        {
                        next CHAR if ($line1[$j] eq $line2[$j]);
                        #Skip positions with both sequences containing gaps
(in case we haven't filtered them before)
                        if ($extend==0)
                                {
                                $score-=$open_penalty;
                                $extend=1;
                                }
                        else {$score-=$extend_penalty;}
                        }
                else
                        {
                        $score+=$matrix->get_entry($line1[$j],$line2[$j]);
                        $extend=0;
                        }
                }
        return $score;
        }


2014-10-01 1:54 GMT+09:00 Fields, Christopher J <cjfields at illinois.edu>:

>  On Sep 30, 2014, at 1:20 AM, Alexey Morozov <alexeymorozov1991 at gmail.com>
> wrote:
>
>  Dear colleagues,
>  Is there a method in bioperl that will calculate pairwise alignment
> scores for any given pair of genes in MSA (according to a given matrix and
> gap opening/extension cost)? It seems that Bio::SimpleAlign methods only
> work with score if it has been described in MSA file and can only hold a
> general multiple sequence alignment score.
>
> --
> Alexey Morozov,
> LIN SB RAS, bioinformatics group.
> Irkutsk, Russia.
>
>
>  Bio::Align::PairwiseStatistics appears to deal with this, see if it fits
> your needs.  You may need to extract the pairwise alignment of the two
> genes if they are in a multiple sequence alignment.
>
>  chris
>
>


-- 
Alexey Morozov,
LIN SB RAS, bioinformatics group.
Irkutsk, Russia.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.open-bio.org/pipermail/bioperl-l/attachments/20141001/3f93fda3/attachment.html>


More information about the Bioperl-l mailing list