Open main menu

Biolecture.org β

MSA code 1

Revision as of 08:22, 7 December 2018 by imported>Na kyung Jung (Created page with "<p>Due to my lack of coding ability, I couldn't solve some errors.</p> <p>#!/usr/bin/perl</p> <p> </p> <p>use 5.010;</p> <p>use strict;</p> <p>use warnings;</p> <p...")
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)

Due to my lack of coding ability, I couldn't solve some errors.

#!/usr/bin/perl

 

use 5.010;

use strict;

use warnings;

 

my $seq1= 'AAGAATAGTATTTCGCTTTTTTATA';

my $seq2= 'AGAAATAGTATTTCGGTTAATTATA';

my $seq3= 'AGACATACTAATTCGGTCCATTATA';

my @matrix;

$matrix[0][0]{score}   = 0;

$matrix[0][0]{pointer} = "none";

for(my $j = 1; $j <= length($seq1); $j++) {

    $matrix[0][$j]{score}   = $GAP * $j;

    $matrix[0][$j]{pointer} = "left";

}

for (my $i = 1; $i <= length($seq2); $i++) {

    $matrix[$i][0]{score}   = $GAP * $i;

    $matrix[$i][0]{pointer} = "up";

}

 

# fill

for(my $i = 1; $i <= length($seq2); $i++) {

    for(my $j = 1; $j <= length($seq1); $j++) {

        my ($diagonal_score, $left_score, $up_score);

 

        # calculate match score

        my $letter1 = substr($seq1, $j-1, 1);

        my $letter2 = substr($seq2, $i-1, 1);                        

+  

        if ($letter1 eq $letter2) {

            $diagonal_score = $matrix[$i-1][$j-1]{score} + $MATCH;

        }

        else {

            $diagonal_score = $matrix[$i-1][$j-1]{score} + $MISMATCH;

        }

 

        # calculate gap scores

        $up_score   = $matrix[$i-1][$j]{score} + $GAP;

        $left_score = $matrix[$i][$j-1]{score} + $GAP;

 

        # choose best score

        if ($diagonal_score >= $up_score) {

            if ($diagonal_score >= $left_score) {

                $matrix[$i][$j]{score}   = $diagonal_score;

                $matrix[$i][$j]{pointer} = "diagonal";

            }

        else {

                $matrix[$i][$j]{score}   = $left_score;

                $matrix[$i][$j]{pointer} = "left";

            }

        } else {

            if ($up_score >= $left_score) {

                $matrix[$i][$j]{score}   = $up_score;

                $matrix[$i][$j]{pointer} = "up";

            }

            else {

                $matrix[$i][$j]{score}   = $left_score;

                $matrix[$i][$j]{pointer} = "left";

            }

        }

    }

}

# trace-back

 

my $align1 = "";

my $align2 = "";

 

# start at last cell of matrix

my $j = length($seq1);

my $i = length($seq2);

 

while (1) {

    last if $matrix[$i][$j]{pointer} eq "none"; # ends at first cell o

+f matrix

 

    if ($matrix[$i][$j]{pointer} eq "diagonal") {

        $align1 .= substr($seq1, $j-1, 1);

        $align2 .= substr($seq2, $i-1, 1);

        $i--;

        $j--;

    }

    elsif ($matrix[$i][$j]{pointer} eq "left") {

        $align1 .= substr($seq1, $j-1, 1);

        $align2 .= "-";

        $j--;

    }

    elsif ($matrix[$i][$j]{pointer} eq "up") {

        $align1 .= "-";

        $align2 .= substr($seq2, $i-1, 1);

        $i--;

    }   

}

$align1 = reverse $align1;

$align2 = reverse $align2;

print "$align1\n";

print "$align2\n";

 

$matrix[0][0]{score}   = 0;

$matrix[0][0]{pointer} = "none";

for(my $j = 1; $j <= length($seq1); $j++) {

    $matrix[0][$j]{score}   = $GAP * $j;

    $matrix[0][$j]{pointer} = "left";

}

for (my $i = 1; $i <= length($seq3); $i++) {

    $matrix[$i][0]{score}   = $GAP * $i;

    $matrix[$i][0]{pointer} = "up";

}

 

# fill

for(my $i = 1; $i <= length($seq3); $i++) {

    for(my $j = 1; $j <= length($seq1); $j++) {

        my ($diagonal_score, $left_score, $up_score);

 

        # calculate match score

        my $letter1 = substr($align1, $j-1, 1);

        my $letter3 = substr($seq3, $i-1, 1);                        

+  

        if ($letter1 eq $letter3) {

            $diagonal_score = $matrix[$i-1][$j-1]{score} + $MATCH;

        }

        else {

            $diagonal_score = $matrix[$i-1][$j-1]{score} + $MISMATCH;

        }

 

        # calculate gap scores

        $up_score   = $matrix[$i-1][$j]{score} + $GAP;

        $left_score = $matrix[$i][$j-1]{score} + $GAP;

 

        # choose best score

        if ($diagonal_score >= $up_score) {

            if ($diagonal_score >= $left_score) {

                $matrix[$i][$j]{score}   = $diagonal_score;

                $matrix[$i][$j]{pointer} = "diagonal";

            }

        else {

                $matrix[$i][$j]{score}   = $left_score;

                $matrix[$i][$j]{pointer} = "left";

            }

        } else {

            if ($up_score >= $left_score) {

                $matrix[$i][$j]{score}   = $up_score;

                $matrix[$i][$j]{pointer} = "up";

            }

            else {

                $matrix[$i][$j]{score}   = $left_score;

                $matrix[$i][$j]{pointer} = "left";

            }

        }

    }

}

# trace-back

 

my $alignalign1 = "";

my $align3 = "";

 

# start at last cell of matrix

my $j = length($seq1);

my $i = length($seq3);

 

while (1) {

    last if $matrix[$i][$j]{pointer} eq "none"; # ends at first cell o

+f matrix

 

    if ($matrix[$i][$j]{pointer} eq "diagonal") {

        $alignalign1 .= substr($align1, $j-1, 1);

        $align3 .= substr($seq3, $i-1, 1);

        $i--;

        $j--;

    }

    elsif ($matrix[$i][$j]{pointer} eq "left") {

        $alignalign1 .= substr($align1, $j-1, 1);

        $align3 .= "-";

        $j--;

    }

    elsif ($matrix[$i][$j]{pointer} eq "up") {

        $alignalign1 .= "-";

        $align3 .= substr($seq3, $i-1, 1);

        $i--;

    }   

}

$alignalign1 = reverse $alignalign1;

$align3 = reverse $align3;

print "$alignalign1\n";

print "$align3\n";

 

$matrix[0][0]{score}   = 0;

$matrix[0][0]{pointer} = "none";

for(my $j = 1; $j <= length($seq1); $j++) {

    $matrix[0][$j]{score}   = $GAP * $j;

    $matrix[0][$j]{pointer} = "left";

}

for (my $i = 1; $i <= length($seq3); $i++) {

    $matrix[$i][0]{score}   = $GAP * $i;

    $matrix[$i][0]{pointer} = "up";

}

 

# fill

for(my $i = 1; $i <= length($seq3); $i++) {

    for(my $j = 1; $j <= length($seq2); $j++) {

        my ($diagonal_score, $left_score, $up_score);

 

        # calculate match score

        my $letter2 = substr($align2, $j-1, 1);

        my $letter3 = substr($seq3, $i-1, 1);                        

+  

        if ($letter2 eq $letter3) {

            $diagonal_score = $matrix[$i-1][$j-1]{score} + $MATCH;

        }

        else {

            $diagonal_score = $matrix[$i-1][$j-1]{score} + $MISMATCH;

        }

 

        # calculate gap scores

        $up_score   = $matrix[$i-1][$j]{score} + $GAP;

        $left_score = $matrix[$i][$j-1]{score} + $GAP;

 

        # choose best score

        if ($diagonal_score >= $up_score) {

            if ($diagonal_score >= $left_score) {

                $matrix[$i][$j]{score}   = $diagonal_score;

                $matrix[$i][$j]{pointer} = "diagonal";

            }

        else {

                $matrix[$i][$j]{score}   = $left_score;

                $matrix[$i][$j]{pointer} = "left";

            }

        } else {

            if ($up_score >= $left_score) {

                $matrix[$i][$j]{score}   = $up_score;

                $matrix[$i][$j]{pointer} = "up";

            }

            else {

                $matrix[$i][$j]{score}   = $left_score;

                $matrix[$i][$j]{pointer} = "left";

            }

        }

    }

}

# trace-back

 

my $alignalign2 = "";

my $align3 = "";

 

# start at last cell of matrix

my $j = length($seq2);

my $i = length($seq3);

 

while (1) {

    last if $matrix[$i][$j]{pointer} eq "none"; # ends at first cell o

+f matrix

 

    if ($matrix[$i][$j]{pointer} eq "diagonal") {

        $alignalign2 .= substr($align2, $j-1, 1);

        $align3 .= substr($seq3, $i-1, 1);

        $i--;

        $j--;

    }

    elsif ($matrix[$i][$j]{pointer} eq "left") {

        $alignalign2 .= substr($align2, $j-1, 1);

        $align3 .= "-";

        $j--;

    }

    elsif ($matrix[$i][$j]{pointer} eq "up") {

        $alignalign2 .= "-";

        $align3 .= substr($seq3, $i-1, 1);

        $i--;

    }   

}

$alignalign2 = reverse $alignalign2;

$align3 = reverse $align3;

print "$alignalign1\n";

print "$alignalign2\n";

print "$align3\n";