5.3 Restriction.pm: Finding Recognition Sites
The
time has now come to write the class
that creates an object out of a sequence, enzyme name(s), and the map
of the location(s) of the enzyme recognition sites in the sequence.
This module depends a great deal on the module
Rebase developed in the previous section, but it
is fairly short because it just tries to do a small job. This new
Restriction class takes a
Rebase object (which has the Rebase database
translated to regular expressions), some sequence, and a list of
enzymes, and uses the regular expressions to find the recognition
sites in the sequence. (Note that it doesn't use
inheritance; it simply creates a Rebase object to
use.)
In this module, the restriction map (the list of locations where the
enzymes have recognition sites in the sequence) is obtained through
an auxiliary method map_enzyme that simply lists
the locations. Clearly, a more graphical display would be easier and
more useful. I'll consider that as the book
progresses.
5.3.1 The Restriction.pm Module
Here is the Restriction.pm module:
package Restriction;
#
# A class to find locations of restriction enzyme recognition sites in
# DNA sequence data.
#
use strict;
use warnings;
use Carp;
# Class data and methods
{
# A list of all attributes with default values.
# "enzyme" is given as an argument possibly multiple time, set as key to _map hash
my %_attributes = (
_rebase => { }, # A Rebase.pm hash-based object
# key = restriction enzyme name
# value = space-separated string of recognition sites => regular expressions
_sequence => '', # DNA sequence data in raw format (only bases)
_map => { },# a hash: keys are enzyme names,
# values are arrays of locations
_enzyme => '', # space- or comma-separated enzyme names,
# set as key to _map hash
);
# Global variable to keep count of existing objects
my $_count = 0;
# Return a list of all attributes
sub _all_attributes {
keys %_attributes;
}
# Manage the count of existing objects
sub get_count {
$_count;
}
sub _incr_count {
++$_count;
}
sub _decr_count {
--$_count;
}
}
This opening block shows that Restriction is also
a class with a small and simple set of attributes. One is just the
DNA sequence data, _sequence.
The attribute _map is a hash that stores the
computed restriction map with a key for each restriction enzyme and a
value consisting of an array of the positions in the sequence where
recognition sites for that enzyme occur.
The attribute _enzyme is one or more enzyme names
separated by spaces or commas.
The remaining attribute _rebase is an object of
the class I developed in the last section, the class
Rebase. Recall that a Rebase
object gives you the ability to get regular expressions for
recognition sites of restriction enzymes. When you call
Restriction->new you must include as one of its
arguments a Rebase object, as is shown in the
example given in the POD documentation later in this section.
Since a program can have many restriction maps, I also maintain the
count of Restriction objects.
5.3.1.1 Initializing Restriction objects
Next, the new
constructor method creates and initializes
Restriction objects:
# The "new" constructor method, called from class, e.g.
sub new {
my ($class, %args) = @_;
# Create a new object
my $self = bless { }, $class;
# Set the attributes for the provided arguments
foreach my $attribute ($self->_all_attributes( )) {
# E.g. attribute = "_name", argument = "name"
my($argument) = ($attribute =~ /^_(.*)/);
if (exists $args{$argument}) {
if($argument eq 'enzyme') {
# permit space or comma separated enzyme names
$args{$argument} =~ s/,/ /g;
}
$self->{$attribute} = $args{$argument};
}
}
# Check that the correct arguments are given
if( not defined $self->{_rebase} ) {
croak "A Rebase object must be given as an argument";
}elsif( ref($self->{_rebase}) ne 'Rebase' ) {
croak "The argument to rebase is not a Rebase object";
}elsif( not defined $self->{_sequence} ) {
croak "A sequence must be given as an argument";
}
# Calculate the locations for each enzyme, store in _map hash attribute
foreach my $enzyme (split(" ", $self->{_enzyme})) {
$self->map_enzyme($enzyme);
}
$self->_incr_count;
return $self;
}
# For this simple class I have no AUTOLOAD or DESTROY
# No get_rebase method, I don't want to pass around a huge hash
# No set mutators: all initialization done by way of "new" constructor
# No clone method. Each sequence and set of enzymes can be easily calculated
# by means of a "new" command.
5.3.1.2 The methods explained
The new constructor method checks that the proper
and required arguments are provided: a sequence, enzyme(s), and a
Rebase object.
Then the new constructor performs the required
computation by calling the method map_enzyme for
each enzyme to determine the (possibly empty) array of locations in
which the enzyme can cut the sequence. These enzyme-array key/value
pairs are stored in the _map attribute.
As you can see in the next segment of code from
Restriction.pm, the method
map_enzyme works by getting the regular
expressions that have been calculated for the enzyme and then finding
the positions where the regular expressions match the sequence.
The get_regular_expressions method accesses the
Rebase object that was passed to the
Restriction->new constructor method and appears
as the $self->{_rebase} attribute. In that
object, the attribute _rebase is the hash of the
database, which looks for the value of the key
$enzyme.
The match_positions method puts a match of a
regular expression in a while loop, where it loops
once for each location it finds. The offset of the matching sequence
is available in the special variable $-[0] after a
successful regular-expression match (see the
perlvar section of the Perl documentation for more
details).
sub map_enzyme {
my($self, $enzyme) = @_;
my(@positions) = ( );
my(@res) = $self->get_regular_expressions($enzyme);
foreach my $re (@res) {
push @positions, $self->match_positions($re);
}
@{$self->{_map}{$enzyme}} = @positions;
return @positions;
}
sub get_regular_expressions {
my($self, $enzyme) = @_;
my(%sites) = split(' ', $self->{_rebase}{_rebase}{$enzyme});
# May have duplicate values
return values %sites;
}
# Find positions of a regular expression in the sequence
sub match_positions {
my($self, $regexp) = @_;
my @positions = ( );
# Determine positions of regular expression matches
while ( $self->{_sequence} =~ /$regexp/ig ) {
push @positions, ($-[0] + 1 );
}
return(@positions);
}
Finally, here are some helper methods that report on the attributes
of sequence and of the calculated map for a given enzyme (or the
entire hash that contains the map of all the enzymes). I
didn't use AUTOLOAD here because
there are only a handful of attributes that return a few different
data types: array, scalar string, and hash.
sub get_enzyme_map {
my($self, $enzyme) = @_;
@{$self->{_map}{$enzyme}};
}
sub get_enzyme_names {
my($self) = @_;
keys %{$self->{_map}};
}
sub get_sequence {
my($self) = @_;
$self->{_sequence};
}
sub get_map {
my($self) = @_;
%{$self->{_map}};
}
5.3.1.3 Documentation
Here's
the (bare bones) documentation for the class embedded right in the
module using the POD documentation language:
=head1 Restriction
Restriction: Given a Rebase object, sequence, and list of restriction enzyme
names, return the locations of the recognition sites in the sequence
=head1 Synopsis
use Restriction;
use Rebase;
use strict;
use warnings;
my $rebase = Rebase->new(
dbmfile => 'BIONET',
bionetfile => 'bionet.212'
);
my $restrict = Restriction->new(
rebase => $rebase,
enzyme => 'EcoRI, HindIII',
sequence => 'ACGAATTCCGGAATTCG',
);
print "Locations for EcoRI are ", join(' ',
$restrict->get_enzyme_map('EcoRI')), "\n";
=head1 AUTHOR
James Tisdall
=head1 COPYRIGHT
Copyright (c) 2003, James Tisdall
=cut
1;
Saving the Synopsis example in a file and running
it (making sure the required bionet.212 file is in
place) gives the output:
EcoRI data in Rebase is GAATTC GAATTC
Sequence is ACGAATTCCGGAATTCG
Locations for EcoRI are 3 11
|
لطفا منتظر باشید ...
