Bio/Tools/Gel.pm

   1 #
   2 # BioPerl module for Bio::Tools::Gel
   3 # Copyright Allen Day <allenday@ucla.edu>
   4 # You may distribute this module under the same terms as perl itself
   5
   6 # POD documentation - main docs before the code
   7
   8 =head1 NAME
   9
  10 Bio::Tools::Gel - Calculates relative electrophoretic migration distances
  11
  12 =head1 SYNOPSIS
  13
  14     use Bio::PrimarySeq;
  15     use Bio::Restriction::Analysis;
  16     use Bio::Tools::Gel;
  17
  18     # get a sequence
  19     my $d = 'AAAAAAAAAGAATTCTTTTTTTTTTTTTTGAATTCGGGGGGGGGGGGGGGGGGGG';
  20     my $seq1 = Bio::Seq->new(-id=>'groundhog day',-seq=>$d);
  21
  22     # cut it with an enzyme
  23     my $ra=Bio::Restriction::Analysis->new(-seq=>$seq1);
  24     @cuts = $ra->fragments('EcoRI'), 3;
  25
  26     # analyse the fragments in a gel
  27     my $gel = Bio::Tools::Gel->new(-seq=>\@cuts,-dilate=>10);
  28     my %bands = $gel->bands;
  29     foreach my $band (sort {$b <=> $a} keys %bands){
  30       print $band,"\t", sprintf("%.1f", $bands{$band}),"\n";
  31     }
  32
  33     #prints:
  34     #20   27.0
  35     #25   26.0
  36     #10   30.0
  37
  38
  39 =head1 DESCRIPTION
  40
  41 This takes a set of sequences or Bio::Seq objects, and calculates their
  42 respective migration distances using:
  43     distance = dilation * (4 - log10(length(dna));
  44
  45 Source: Molecular Cloning, a Laboratory Manual. Sambrook, Fritsch, Maniatis.
  46 CSHL Press, 1989.
  47
  48 Bio::Tools::Gel currently calculates migration distances based solely on
  49 the length of the nucleotide sequence.  Secondary or tertiary structure,
  50 curvature, and other biophysical attributes of a sequence are currently
  51 not considered.  Polypeptide migration is currently not supported.
  52
  53 =head1 FEEDBACK
  54
  55 =head2 Mailing Lists
  56
  57 User feedback is an integral part of the evolution of this and other
  58 Bioperl modules. Send your comments and suggestions preferably to
  59 the Bioperl mailing list.  Your participation is much appreciated.
  60
  61   bioperl-l@bioperl.org                  - General discussion
  62   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  63
  64 =head2 Support
  65
  66 Please direct usage questions or support issues to the mailing list:
  67
  68 I<bioperl-l@bioperl.org>
  69
  70 rather than to the module maintainer directly. Many experienced and
  71 reponsive experts will be able look at the problem and quickly
  72 address it. Please include a thorough description of the problem
  73 with code and data examples if at all possible.
  74
  75 =head2 Reporting Bugs
  76
  77 Report bugs to the Bioperl bug tracking system to help us keep track
  78 of the bugs and their resolution. Bug reports can be submitted via the
  79 web:
  80
  81   https://redmine.open-bio.org/projects/bioperl/
  82
  83 =head1 AUTHOR - Allen Day
  84
  85 Email allenday@ucla.edu
  86
  87 =head1 APPENDIX
  88
  89 The rest of the documentation details each of the object methods.
  90 Internal methods are usually preceded with a _
  91
  92 =cut
  93
  94
  95 # Let the code begin...
  96
  97
  98 package Bio::Tools::Gel;
  99 use strict;
 100
 101 use Bio::PrimarySeq;
 102
 103 use base qw(Bio::Root::Root);
 104
 105 =head2 new
 106
 107  Title   : new
 108  Usage   : my $gel = Bio::Tools::Gel->new(-seq => $sequence,-dilate => 3);
 109  Function: Initializes a new Gel
 110  Returns : Bio::Tools::Gel
 111  Args    : -seq      => Bio::Seq(s), scalar(s) or list of either/both
 112                         (default: none)
 113            -dilate   => Expand band migration distances (default: 1)
 114
 115 =cut
 116
 117 sub new {
 118   my($class,@args) = @_;
 119
 120   my $self = $class->SUPER::new(@args);
 121   my ($seqs,$dilate) = $self->_rearrange([qw(SEQ DILATE)],
 122                                           @args);
 123   if( ! ref($seqs)  ) {
 124       $self->add_band([$seqs]);
 125   } elsif( ref($seqs) =~ /array/i ||
 126            $seqs->isa('Bio::PrimarySeqI') ) {
 127       $self->add_band($seqs);
 128   }
 129   $self->dilate($dilate || 1);
 130
 131   return $self;
 132 }
 133
 134
 135 =head2 add_band
 136
 137  Title   : add_band
 138  Usage   : $gel->add_band($seq);
 139  Function: Calls _add_band with a (possibly created) Bio::Seq object.
 140  Returns :
 141  Args    : Bio::Seq, scalar sequence, or list of either/both.
 142
 143 =cut
 144
 145 sub add_band {
 146   my($self,$args) = @_;
 147
 148   foreach my $arg (@$args){
 149       my $seq;
 150       if( ! ref($arg) ) {
 151           if( $arg =~ /^\d+/ ) {
 152               $seq= Bio::PrimarySeq->new(-seq=>"N"x$arg, -id => $arg);
 153           } else {
 154               $seq= Bio::PrimarySeq->new(-seq=>$arg,-id=>length($arg));
 155           }
 156       } elsif( $arg->isa('Bio::PrimarySeqI') ) {
 157           $seq = $arg;
 158       }
 159
 160     $seq->validate_seq or $seq->throw("invalid symbol in sequence".$seq->seq()."\n");
 161     $self->_add_band($seq);
 162   }
 163 }
 164
 165 =head2 _add_band
 166
 167  Title   : _add_band
 168  Usage   : $gel->_add_band($seq);
 169  Function: Adds a new band to the gel.
 170  Returns :
 171  Args    : Bio::Seq object
 172
 173 =cut
 174
 175 sub _add_band {
 176   my($self,$arg) = @_;
 177   if( defined $arg) {
 178       push (@{$self->{'bands'}},$arg);
 179   }
 180 }
 181
 182 =head2 dilate
 183
 184  Title   : dilate
 185  Usage   : $gel->dilate(1);
 186  Function: Sets/retrieves the dilation factor.
 187  Returns : dilation factor
 188  Args    : Float or none
 189
 190 =cut
 191
 192 sub dilate {
 193   my($self,$arg) = @_;
 194   return $self->{dilate} unless $arg;
 195   $self->throw("-dilate should be numeric") if defined $arg and $arg =~ /[^e\d\.]/;
 196   $self->{dilate} = $arg;
 197   return $self->{dilate};
 198 }
 199
 200 sub migrate {
 201   my ($self,$arg) = @_;
 202   $arg = $self unless $arg;
 203   if ( $arg ) {
 204       return 4 - log10($arg);
 205   } else { return 0; }
 206 }
 207
 208 =head2 bands
 209
 210  Title   : bands
 211  Usage   : $gel->bands;
 212  Function: Calculates migration distances of sequences.
 213  Returns : hash of (seq_id => distance)
 214  Args    :
 215
 216 =cut
 217
 218 sub bands {
 219   my $self = shift;
 220   $self->throw("bands() is read-only") if @_;
 221
 222   my %bands = ();
 223
 224   foreach my $band (@{$self->{bands}}){
 225       my $distance = $self->dilate * migrate($band->length);
 226       $bands{$band->id} = $distance;
 227   }
 228
 229   return %bands;
 230 }
 231
 232 =head2 log10
 233
 234  Title   : log10
 235  Usage   : log10($n);
 236  Function: returns base 10 log of $n.
 237  Returns : float
 238  Args    : float
 239
 240 =cut
 241
 242 #from programming perl
 243 sub log10 {
 244     my $n = shift;
 245     return log($n)/log(10);
 246 }
 247
 248 1;