Re: Getting pairwise alignment scores for existing multiple alignment

Alexey Morozov <[email protected]>
Newsgroups gmane.comp.lang.perl.bio.general
Message-ID <CAK=7cMMwgUwoAsLjY5gTnNMNKMKCMqUtNnc1Pauja3eHCrBkHg@mail.gmail.com>
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 <[email protected]>:

>  On Sep 30, 2014, at 1:20 AM, Alexey Morozov <[email protected]>
> 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.

_______________________________________________
Bioperl-l mailing list
[email protected]
http://mailman.open-bio.org/mailman/listinfo/bioperl-l
lmpx.com only provides a reader for public news (NNTP) servers. It is not affiliated with the servers or forums shown here and is not responsible for the content of articles, which is written by their respective authors.