~ubuntu-branches/ubuntu/hoary/bioperl/hoary

« back to all changes in this revision

Viewing changes to Bio/Seq/Meta/Array.pm

  • Committer: Bazaar Package Importer
  • Author(s): Matt Hope
  • Date: 2004-04-18 14:24:11 UTC
  • mfrom: (1.2.1 upstream) (2.1.1 warty)
  • Revision ID: james.westby@ubuntu.com-20040418142411-gr92uexquw4w8liq
Tags: 1.4-1
* New upstream release
* Examples and working code are installed by default to usr/bin,
  this has been moved to usr/share/doc/bioperl/bin

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
# $Id: Array.pm,v 1.5 2003/09/08 12:17:14 heikki Exp $
 
2
#
 
3
# BioPerl module for Bio::Seq::Meta::Array
 
4
#
 
5
# Cared for by Heikki Lehvaslaiho
 
6
#
 
7
# Copyright Heikki Lehvaslaiho
 
8
#
 
9
# You may distribute this module under the same terms as perl itself
 
10
 
 
11
# POD documentation - main docs before the code
 
12
 
 
13
=head1 NAME
 
14
 
 
15
Bio::Seq::Meta::Array - array-based generic implementation of a sequence class with residue-based meta information
 
16
 
 
17
=head1 SYNOPSIS
 
18
 
 
19
  use Bio::LocatableSeq;
 
20
  use Bio::Seq::Meta::Array;
 
21
 
 
22
  my $seq = Bio::LocatableSeq->new(-id=>'test',
 
23
                                   -seq=>'ACTGCTAGCT',
 
24
                                   -start=>2434,
 
25
                                   -start=>2443,
 
26
                                   -strand=>1,
 
27
                                   -varbose=>1, # to see warnings
 
28
                                  );
 
29
  bless $seq, Bio::Seq::Meta::Array;
 
30
  # the existing sequence object can be a Bio::PrimarySeq, too
 
31
 
 
32
  # to test this is a meta seq object
 
33
  $seq->isa("Bio::Seq::Meta::Array")
 
34
      || $seq->throw("$seq is not a Bio::Seq::Meta::Array");
 
35
 
 
36
  $seq->meta('1 2 3 4 5 6 7 8 9 10');
 
37
 
 
38
  # accessors
 
39
  $arrayref   = $seq->meta();
 
40
  $string     = $seq->meta_text();
 
41
  $substring  = $seq->submeta_text(2,5);
 
42
  $unique_key = $seq->accession_number();
 
43
 
 
44
=head1 DESCRIPTION
 
45
 
 
46
This class implements generic methods for sequences with residue-based
 
47
meta information. Meta sequences with meta data are Bio::LocatableSeq
 
48
objects with additional methods to store that meta information. See
 
49
L<Bio::LocatableSeq> and L<Bio::Seq::MetaI>.
 
50
 
 
51
The meta information in this class can be a string of variable length
 
52
and can be a complex structure.  Blank values are undef or zero.
 
53
 
 
54
Application specific implementations should inherit from this class to
 
55
override and add to these methods.
 
56
 
 
57
=head1 SEE ALSO
 
58
 
 
59
L<Bio::LocatableSeq>, 
 
60
L<Bio::Seq::MetaI>, 
 
61
L<Bio::Seq::Meta>, 
 
62
 
 
63
=head1 FEEDBACK
 
64
 
 
65
=head2 Mailing Lists
 
66
 
 
67
User feedback is an integral part of the evolution of this and other
 
68
Bioperl modules. Send your comments and suggestions preferably to one
 
69
of the Bioperl mailing lists.  Your participation is much appreciated.
 
70
 
 
71
  bioperl-l@bioperl.org                      - General discussion
 
72
  http://bio.perl.org/MailList.html          - About the mailing lists
 
73
 
 
74
=head2 Reporting Bugs
 
75
 
 
76
Report bugs to the Bioperl bug tracking system to help us keep track
 
77
the bugs and their resolution.  Bug reports can be submitted via email
 
78
or the web:
 
79
 
 
80
  bioperl-bugs@bio.perl.org
 
81
  http://bugzilla.bioperl.org/
 
82
 
 
83
=head1 AUTHOR - Heikki Lehvaslaiho
 
84
 
 
85
Email heikki@ebi.ac.uk
 
86
 
 
87
=head1 CONTRIBUTORS
 
88
 
 
89
Chad Matsalla, bioinformatics@dieselwurks.com
 
90
Aaron Mackey, amackey@virginia.edu
 
91
 
 
92
=head1 APPENDIX
 
93
 
 
94
The rest of the documentation details each of the object methods.
 
95
Internal methods are usually preceded with a _
 
96
 
 
97
=cut
 
98
 
 
99
 
 
100
# Let the code begin...
 
101
 
 
102
 
 
103
package Bio::Seq::Meta::Array;
 
104
use vars qw(@ISA $DEFAULT_NAME $GAP $META_GAP);
 
105
use strict;
 
106
use Bio::LocatableSeq;
 
107
use Bio::Seq::MetaI;
 
108
 
 
109
use Data::Dumper;
 
110
 
 
111
#use overload '""' => \&to_string;
 
112
 
 
113
@ISA = qw( Bio::LocatableSeq Bio::Seq Bio::Seq::MetaI );
 
114
 
 
115
BEGIN {
 
116
 
 
117
    $DEFAULT_NAME = 'DEFAULT';
 
118
    $GAP = '-';
 
119
    $META_GAP = 0;
 
120
}
 
121
 
 
122
=head2 new
 
123
 
 
124
 Title   : new
 
125
 Usage   : $metaseq = Bio::Seq::Meta->new
 
126
                ( -meta => 'aaaaaaaabbbbbbbb',
 
127
                  -seq =>  'TKLMILVSHIVILSRM'
 
128
                  -id  => 'human_id',
 
129
                  -accession_number => 'S000012',
 
130
                );
 
131
 Function: Constructor for Bio::Seq::Meta class, meta data being in a
 
132
           string. Note that you can provide an empty quality string.
 
133
 Returns : a new Bio::Seq::Meta object
 
134
 
 
135
=cut
 
136
 
 
137
sub new {
 
138
    my ($class, %args) = @_;
 
139
        #defined inheritance according to stated baseclass,
 
140
        #if undefined then will be PrimarySeq
 
141
        if (defined($args{'-baseclass'})) {
 
142
                @ISA = ($args{'-baseclass'},"Bio::Seq::MetaI");
 
143
                }
 
144
        else {
 
145
                @ISA = qw( Bio::LocatableSeq Bio::Seq Bio::Seq::MetaI );
 
146
                }
 
147
 
 
148
    my $self = $class->SUPER::new(%args);
 
149
 
 
150
    my($meta) =
 
151
        $self->_rearrange([qw(META
 
152
                              )],
 
153
                          %args);
 
154
 
 
155
    $self->{'_meta'}->{$DEFAULT_NAME} = undef;
 
156
 
 
157
    $meta && $self->meta($meta);
 
158
 
 
159
    return $self;
 
160
}
 
161
 
 
162
 
 
163
=head2 meta
 
164
 
 
165
 Title   : meta
 
166
 Usage   : $meta_values  = $obj->meta($values_string);
 
167
 Function:
 
168
 
 
169
           Get and set method for the meta data starting from residue
 
170
           position one. Since it is dependent on the length of the
 
171
           sequence, it needs to be manipulated after the sequence.
 
172
 
 
173
           The length of the returned value always matches the length
 
174
           of the sequence.
 
175
 
 
176
 Returns : reference to an array of meta data
 
177
 Args    : new value, string or array ref, optional
 
178
 
 
179
=cut
 
180
 
 
181
sub meta {
 
182
   shift->named_meta($DEFAULT_NAME, shift);
 
183
}
 
184
 
 
185
=head2 meta_text
 
186
 
 
187
 Title   : meta_text
 
188
 Usage   : $meta_values  = $obj->meta_text($values_arrayref);
 
189
 Function: Variant of meta() guarantied to return a textual
 
190
           representation  of meta data. For details, see L<meta>.
 
191
 Returns : a string
 
192
 Args    : new value, string or array ref, optional
 
193
 
 
194
=cut
 
195
 
 
196
sub meta_text {
 
197
    return join ' ',  @{shift->meta(shift)};
 
198
}
 
199
 
 
200
=head2 named_meta
 
201
 
 
202
 Title   : named_meta()
 
203
 Usage   : $meta_values  = $obj->named_meta($name, $values_arrayref);
 
204
 Function: A more general version of meta(). Each meta data set needs
 
205
           to be named. See also L<meta_names>.
 
206
 Returns : reference to an array of meta data
 
207
 Args    : scalar, name of the meta data set
 
208
           new value, string or array ref, optional
 
209
 
 
210
=cut
 
211
 
 
212
sub named_meta {
 
213
   my ($self, $name, $value) = @_;
 
214
 
 
215
   $name ||= $DEFAULT_NAME;
 
216
 
 
217
   if (defined $value) {
 
218
       my ($arrayref);
 
219
 
 
220
       if (ref $value eq 'ARRAY' ) { # array ref
 
221
           $arrayref = $value;
 
222
       }
 
223
       elsif (not ref($value)) { # scalar
 
224
           $arrayref = [split /\s+/, $value];
 
225
       } else {
 
226
           $self->throw("I need a scalar or array ref, not [". ref($value). "]");
 
227
       }
 
228
 
 
229
       # test for length
 
230
       my $diff = $self->length - @{$arrayref};
 
231
       if ($diff > 0) {
 
232
           foreach (1..$diff) { push @{$arrayref}, 0;}
 
233
       }
 
234
 
 
235
       $self->{'_meta'}->{$name} = $arrayref;
 
236
 
 
237
       #$self->_test_gap_positions($name) if $self->verbose > 0;
 
238
   }
 
239
 
 
240
   if (defined $self->{'_meta'}->{$name} and
 
241
       scalar @{$self->{'_meta'}->{$name}} > $self->length ) {
 
242
       return [@{$self->{'_meta'}->{$name}}[0..($self->length-1)]];
 
243
   } else {
 
244
       return $self->{'_meta'}->{$name} || (" " x $self->length);
 
245
   }
 
246
}
 
247
 
 
248
=head2 _test_gap_positions
 
249
 
 
250
 Title   : _test_gap_positions
 
251
 Usage   : $meta_values  = $obj->_test_gap_positions($name);
 
252
 Function: Internal test for correct position of gap characters.
 
253
           Gap being only '-' this time.
 
254
 
 
255
           This method is called from named_meta() when setting meta
 
256
           data but only if verbose is positive as this can be an
 
257
           expensive process on very long sequences. Set verbose(1) to
 
258
           see warnings when gaps do not align in sequence and meta
 
259
           data and turn them into errors by setting verbose(2).
 
260
 
 
261
 Returns : true on success, prints warnings
 
262
 Args    : none
 
263
 
 
264
=cut
 
265
 
 
266
sub _test_gap_positions {
 
267
    my $self = shift;
 
268
    my $name = shift;
 
269
    my $success = 1;
 
270
 
 
271
    $self->seq || return $success;
 
272
    my $len = CORE::length($self->seq);
 
273
    for (my $i=0; $i < $len; $i++) {
 
274
        my $s = substr $self->{seq}, $i, 1;
 
275
        my $m = substr $self->{_meta}->{$name}, $i, 1;
 
276
        $self->warn("Gap mismatch in column [". ($i+1). "] of [$name] meta data in seq [". $self->id. "]")
 
277
            and $success = 0
 
278
                #if ($s eq '-' || $m eq '-') && $s ne $m;
 
279
                if ($m eq '-') && $s ne $m;
 
280
    }
 
281
    return $success;
 
282
}
 
283
 
 
284
=head2 named_meta_text
 
285
 
 
286
 Title   : named_meta_text()
 
287
 Usage   : $meta_values  = $obj->named_meta_text($name, $values_arrayref);
 
288
 Function: Variant of named_meta() guarantied to return a textual
 
289
           representation  of the named meta data.
 
290
           For details, see L<meta>.
 
291
 Returns : a string
 
292
 Args    : scalar, name of the meta data set
 
293
           new value, string or array ref, optional
 
294
 
 
295
=cut
 
296
 
 
297
sub named_meta_text {
 
298
    return join ' ', @{shift->named_meta(@_)};
 
299
 
 
300
}
 
301
 
 
302
=head2 submeta
 
303
 
 
304
 Title   : submeta
 
305
 Usage   : $subset_of_meta_values = $obj->submeta(10, 20, $value_string);
 
306
           $subset_of_meta_values = $obj->submeta(10, undef, $value_string);
 
307
 Function:
 
308
 
 
309
           Get and set method for meta data for subsequences.
 
310
 
 
311
           Numbering starts from 1 and the number is inclusive, ie 1-2
 
312
           are the first two residue of the sequence. Start cannot be
 
313
           larger than end but can be equal.
 
314
 
 
315
           If the second argument is missing the returned values
 
316
           should extend to the end of the sequence.
 
317
 
 
318
           The return value may be a string or an array reference,
 
319
           depending on the implentation. If in doubt, use
 
320
           submeta_text() which is a variant guarantied to return a
 
321
           string.  See L<submeta_text>.
 
322
 
 
323
 Returns : A reference to an array or a string
 
324
 Args    : integer, start position
 
325
           integer, end position, optional when a third argument present
 
326
           new value, string or array ref, optional
 
327
 
 
328
=cut
 
329
 
 
330
sub submeta {
 
331
   shift->named_submeta($DEFAULT_NAME, @_);
 
332
}
 
333
 
 
334
=head2 submeta_text
 
335
 
 
336
 Title   : submeta_text
 
337
 Usage   : $meta_values  = $obj->submeta_text(20, $value_string);
 
338
 Function: Variant of submeta() guarantied to return a textual
 
339
           representation  of meta data. For details, see L<meta>.
 
340
 Returns : a string
 
341
 Args    : new value, string or array ref, optional
 
342
 
 
343
 
 
344
=cut
 
345
 
 
346
sub submeta_text {
 
347
    return join ' ', @{shift->submeta(@_)};
 
348
}
 
349
 
 
350
=head2 named_submeta
 
351
 
 
352
 Title   : named_submeta
 
353
 Usage   : $subset_of_meta_values = $obj->named_submeta($name, 10, 20, $value_string);
 
354
           $subset_of_meta_values = $obj->named_submeta($name, 10);
 
355
 Function: Variant of submeta() guarantied to return a textual
 
356
           representation  of meta data. For details, see L<meta>.
 
357
 Returns : A reference to an array or a string
 
358
 Args    : scalar, name of the meta data set
 
359
           integer, start position
 
360
           integer, end position, optional when a third argument present
 
361
           new value, string or array ref, optional
 
362
 
 
363
=cut
 
364
 
 
365
 
 
366
sub named_submeta {
 
367
    my ($self, $name, $start, $end, $value) = @_;
 
368
 
 
369
    $name ||= $DEFAULT_NAME;
 
370
    $start ||=1;
 
371
    $start =~ /^[+]?\d+$/ and $start > 0 or
 
372
        $self->throw("Need at least a positive integer start value");
 
373
    $start--;
 
374
 
 
375
    my $metaref = $self->{_meta}->{$name};
 
376
 
 
377
    if (defined $value) {
 
378
        my $arrayref;
 
379
 
 
380
        $self->warn("You are setting meta values beyond the length of the sequence\n".
 
381
                    "[$start > ". length($self->seq)."] in sequence ". $self->id)
 
382
            if $start > $self->length;
 
383
 
 
384
        if (ref $value eq 'ARRAY' ) { # array ref
 
385
            $arrayref = $value;
 
386
        }
 
387
        elsif (not ref($value)) { # scalar
 
388
            $arrayref = [split /\s+/, $value];
 
389
        } else {
 
390
            $self->throw("I need a scalar or array ref, not [". ref($value). "]");
 
391
        }
 
392
 
 
393
 
 
394
 
 
395
        $end or $end = @{$arrayref}  + $start;
 
396
        #$end = length @{$arrayref} + $start if $end > (length (@{$arrayref}) + $start);
 
397
        $end--;
 
398
 
 
399
        # test for length; pad if needed
 
400
        my $diff = $end - $start - scalar @{$arrayref};
 
401
        if ($diff > 0) {
 
402
            foreach (1..$diff) { push @{$arrayref}, $META_GAP}
 
403
        }
 
404
 
 
405
        @{$metaref}[$start..$end] = @{$arrayref};
 
406
        return $arrayref;
 
407
 
 
408
    } else {
 
409
 
 
410
        $end or $end = $self->length;
 
411
        $end = $self->length if $end > $self->length;
 
412
        $end--;
 
413
        return [@{$metaref}[$start..$end]];
 
414
 
 
415
    }
 
416
}
 
417
 
 
418
 
 
419
=head2 named_submeta_text
 
420
 
 
421
 Title   : named_submeta_text
 
422
 Usage   : $meta_values  = $obj->named_submeta_text($name, 20, $value_string);
 
423
 Function: Variant of submeta() guarantied to return a textual
 
424
           representation  of meta data. For details, see L<meta>.
 
425
 Returns : a string
 
426
 Args    : scalar, name of the meta data
 
427
 Args    : integer, start position, optional
 
428
           integer, end position, optional
 
429
           new value, string or array ref, optional
 
430
 
 
431
=cut
 
432
 
 
433
sub named_submeta_text {
 
434
    return join ' ', @{shift->named_submeta(@_)};
 
435
}
 
436
 
 
437
=head2 meta_names
 
438
 
 
439
 Title   : meta_names
 
440
 Usage   : @meta_names  = $obj->meta_names()
 
441
 Function: Retrives an array of meta data set names. The default
 
442
           (unnamed) set name is guarantied to be the first name if it
 
443
           contains any data.
 
444
 Returns : an array of names
 
445
 Args    : none
 
446
 
 
447
=cut
 
448
 
 
449
sub meta_names {
 
450
    my ($self) = @_;
 
451
 
 
452
    my @r;
 
453
    foreach  ( sort keys %{$self->{'_meta'}} ) {
 
454
        push (@r, $_) unless $_ eq $DEFAULT_NAME;
 
455
    }
 
456
    unshift @r, $DEFAULT_NAME if $self->{'_meta'}->{$DEFAULT_NAME};
 
457
    return @r;
 
458
 }
 
459
 
 
460
 
 
461
=head1 Bio::PrimarySeqI methods
 
462
 
 
463
=head2 revcom
 
464
 
 
465
 Title   : revcom
 
466
 Usage   : $newseq = $seq->revcom();
 
467
 Function: Produces a new Bio::Seq::MetaI implementing object where
 
468
           the order of residues and their meta information is reversed.
 
469
 Returns : A new (fresh) Bio::Seq::Meta object
 
470
 Args    : none
 
471
 
 
472
=cut
 
473
 
 
474
sub revcom {
 
475
    my $self = shift;
 
476
 
 
477
    my $new = $self->SUPER::revcom;
 
478
    my $end = $self->length - 1;
 
479
    map {
 
480
        $new->{_meta}->{$_} = [ reverse @{$self->{_meta}->{$_}}[0..$end]]
 
481
    } keys %{$self->{_meta}};
 
482
 
 
483
    return $new;
 
484
}
 
485
 
 
486
=head2 trunc
 
487
 
 
488
 Title   : trunc
 
489
 Usage   : $subseq = $seq->trunc(10,100);
 
490
 Function: Provides a truncation of a sequence together with meta data
 
491
 Returns : a fresh Bio::Seq::Meta implementing object
 
492
 Args    : Two integers denoting first and last residue of the sub-sequence.
 
493
 
 
494
=cut
 
495
 
 
496
sub trunc {
 
497
    my ($self, $start, $end) = @_;
 
498
 
 
499
    # test arguments
 
500
    $start =~ /^[+]?\d+$/ and $start > 0 or
 
501
        $self->throw("Need at least a positive integer start value as start");
 
502
    $end =~ /^[+]?\d+$/ and $end > 0 or
 
503
        $self->throw("Need at least a positive integer start value as end");
 
504
    $end >= $start or
 
505
        $self->throw("End position has to be larger or equal to start");
 
506
    $end <= $self->length or
 
507
        $self->throw("End position can not be larger than sequence length");
 
508
 
 
509
 
 
510
    my $new = $self->SUPER::trunc($start, $end);
 
511
    $start--;
 
512
    map {
 
513
        $new->{_meta}->{$_} = [@{$self->{_meta}->{$_}}[$start..$end]]
 
514
    } keys %{$self->{_meta}};
 
515
    return $new;
 
516
}
 
517
 
 
518
1;