Getting pairwise alignment scores for existing multiple alignment

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

Getting pairwise alignment scores for existing multiple alignment

Alexey Morozov-2
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.

_______________________________________________
Bioperl-l mailing list
[hidden email]
http://mailman.open-bio.org/mailman/listinfo/bioperl-l
Reply | Threaded
Open this post in threaded view
|

Re: Getting pairwise alignment scores for existing multiple alignment

Fields, Christopher J
On Sep 30, 2014, at 1:20 AM, Alexey Morozov <[hidden email]> 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


_______________________________________________
Bioperl-l mailing list
[hidden email]
http://mailman.open-bio.org/mailman/listinfo/bioperl-l
Reply | Threaded
Open this post in threaded view
|

Re: Getting pairwise alignment scores for existing multiple alignment

Alexey Morozov-2
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 <[hidden email]>:
On Sep 30, 2014, at 1:20 AM, Alexey Morozov <[hidden email]> 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
[hidden email]
http://mailman.open-bio.org/mailman/listinfo/bioperl-l