如何将其更改为“惯用”Perl?

我也有一些其他的评论,但这里是第一个观察:

my $num_of_rows = length($self->{seq1}) + 1;
my $num_of_columns = length($self->{seq2}) + 1;

所以 $self->{seq1} 和 $self->{seq2} 是字符串,您可以使用 substr . 我希望将它们存储为字符数组:

$self->{seq1} = [ split //, $seq1 ];

以下是我的写作方式:

sub create_matrix {
    my $self = shift;

    my $matrix      = $self->{score_matrix};
    my $path_matrix = $self->{path_matrix};

    my $rows = @{ $self->{seq1} };
    my $cols = @{ $self->{seq2} };

    for my $row (0 .. $rows) {
        $matrix->[$row]->[0] =  $row * $self->{gap_cost};
        $path_matrix->[$row]->[0] = 1;
    }

    my $gap_cost = $self->{gap_cost};

    $matrix->[0] = [ map { $_ * $gap_cost } 0 .. $cols ];
    $path_matrix->[0] = [ (-1) x ($cols + 1) ];

    $path_matrix->[0]->[0] = 2;

    for my $row (1 .. $rows) {
        for my $col (1 .. $cols) {
            my $gap1 = $matrix->[$row - 1]->[$col] + $gap_cost;
            my $gap2 = $matrix->[$row]->[$col - 1] + $gap_cost;
            my $match_mismatch =
                $matrix->[$row - 1]->[$col - 1] +
                $self->get_match_score(
                    $self->{seq1}->[$row - 1],
                    $self->{seq2}->[$col - 1]
                );

            my $max = $matrix->[$row]->[$col] =
                max($gap1, $gap2, $match_mismatch);

            $path_matrix->[$row]->[$col] = $max == $gap1
                    ? -1
                    : $max == $gap2
                    ? 1
                    : 0;
            }
        }
    }

您得到了一些关于语法的建议,但是我也建议使用更模块化的方法,如果没有其他原因的话,代码的可读性。如果你能在担心低级细节之前了解全局,那么就更容易掌握代码的速度。

您的主要方法可能如下所示。

sub create_matrix {
    my $self = shift;
    $self->create_2d_array_of_scores;
    $self->fill_out_first_row;
    $self->fill_out_other_rows;
}

您还可以使用以下几种较小的方法:

n_of_rows
n_of_cols
create_2d_array_of_scores
fill_out_first_row
fill_out_other_rows

通过定义更小的方法(getter、setter等等),您可能会做得更进一步。在那一点上,你的中层方法 create_2d_array_of_scores 根本不会直接接触底层数据结构。

sub matrix      { shift->{score_matrix} }
sub gap_cost    { shift->{gap_cost}     }

sub set_matrix_value {
    my ($self, $r, $c, $val) = @_;
    $self->matrix->[$r][$c] = $val;
}

# Etc.

一个简单的改变是使用 for 这样的循环:

for my $i (0 .. $num_of_rows){
    # Do stuff.
}

有关详细信息,请参阅上的Perl文档 foreach loops 以及 range operator .

而不是像这样取消对二维数组的引用:

$$path_matrix[0][0] = 2;

这样做:

$path_matrix->[0][0] = 2;

此外,您还做了很多if/then/else语句来匹配特定的子序列:这可以更好地写成 given 语句(perl5.10等价于c's switch )请阅读 perldoc perlsyn :

given ($matrix->[$row][$column])
{
    when ($seq1_gap)       { $path_matrix->[$row][$column] = -1; }
    when ($match_mismatch) { $path_matrix->[$row][$column] = 0; }
    when ($seq2_gap)       { $path_matrix->[$row][$column] = 1; }
}

您的大多数代码都在操作二维数组。我认为最大的改进是改用 PDL 如果你想对数组做很多事情,特别是在效率问题上。它是一个Perl模块,提供了出色的数组支持。为了提高效率,底层例程在C中实现,因此速度也很快。

我总是建议你看看 CPAN 对于以前的解决方案或如何在Perl中进行操作的示例。你看过吗 Algorithm::NeedlemanWunsch ?

本模块的文档包括一个匹配DNA序列的示例。下面是一个使用相似矩阵的示例 wikipedia .

#!/usr/bin/perl -w
use strict;
use warnings;
use Inline::Files;                 #multiple virtual files inside code
use Algorithm::NeedlemanWunsch;    # refer CPAN - good style guide

# Read DNA sequences
my @a = read_DNA_seq("DNA_SEQ_A");
my @b = read_DNA_seq("DNA_SEQ_B");

# Read Similarity Matrix (held as a Hash of Hashes)
my %SM = read_Sim_Matrix();

# Define scoring based on "Similarity Matrix" %SM
sub score_sub {
    if ( !@_ ) {
        return -3;                 # gap penalty same as wikipedia)
    }
    return $SM{ $_[0] }{ $_[1] };    # Similarity Value matrix
}

my $matcher = Algorithm::NeedlemanWunsch->new( \&score_sub, -3 );
my $score = $matcher->align( \@a, \@b, { align => \&check_align, } );

print "\nThe maximum score is $score\n";

sub check_align {
    my ( $i, $j ) = @_;              # @a[i], @b[j]
    print "seqA pos: $i, seqB pos: $j\t base \'$a[$i]\'\n";
}

sub read_DNA_seq {
    my $source = shift;
    my @data;
    while (<$source>) {
        push @data, /[ACGT-]{1}/g;
    }
    return @data;
}

sub read_Sim_Matrix {

    #Read DNA similarity matrix (scores per Wikipedia)
    my ( @AoA, %HoH );
    while (<SIMILARITY_MATRIX>) {
        push @AoA, [/(\S+)+/g];
    }

    for ( my $row = 1 ; $row < 5 ; $row++ ) {
        for ( my $col = 1 ; $col < 5 ; $col++ ) {
            $HoH{ $AoA[0][$col] }{ $AoA[$row][0] } = $AoA[$row][$col];
        }
    }
    return %HoH;
}

__DNA_SEQ_A__
A T G T A G T G T A T A G T
A C A T G C A
__DNA_SEQ_B__
A T G T A G T A C A T G C A
__SIMILARITY_MATRIX__
-  A  G  C  T
A  10  -1  -3  -4
G  -1  7  -5  -3
C  -3  -5  9  0
T  -4  -3  0  8

下面是一些示例输出:

seqA pos: 7, seqB pos: 2  base 'G'
seqA pos: 6, seqB pos: 1  base 'T'
seqA pos: 4, seqB pos: 0  base 'A'

The maximum score is 100