BaseFeatureAdaptor.pm 19.3 KB
Newer Older
1
#
2
# EnsEMBL module for Bio::EnsEMBL::DBSQL::BaseAlignFeatureAdaptor
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
#
# Cared for by Ewan Birney <birney@ebi.ac.uk>
#
# Copyright Ewan Birney
#
# You may distribute this module under the same terms as perl itself

# POD documentation - main docs before the code

=head1 NAME

Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor - Abstract Base class for 
                                          FeatureAdaptors

=head1 SYNOPSIS

Abstract class should not be instantiated.  Implementation of
abstract methods must be performed by subclasses.

=head1 DESCRIPTION

This is a base adaptor for feature adaptors. This base class is simply a way
of eliminating code duplication through the implementation of methods 
common to all feature adaptors.

=head1 AUTHOR - Ewan Birney

Email birney@ebi.ac.uk

Describe contact details here

=head1 APPENDIX

The rest of the documentation details each of the object methods. Internal methods are usually preceded with a _

=cut


# Let the code begin...


package Bio::EnsEMBL::DBSQL::BaseFeatureAdaptor;
use vars qw(@ISA);
use strict;

# Object preamble - inherits from Bio::EnsEMBL::Root

use Bio::EnsEMBL::DBSQL::BaseAdaptor;

@ISA = qw(Bio::EnsEMBL::DBSQL::BaseAdaptor);

54
55
56
57
58
59
60
61
62
63
64
65
66
67
68

=head2 generic_fetch

  Arg [1]    : (optional) string $constraint
               An SQL query constraint (i.e. part of the WHERE clause)
  Arg [2]    : (optional) string $logic_name
               the logic_name of the analysis of the features to obtain
  Example    : @fts = $a->generic_fetch('contig_id in (1234, 1235)', 'swall');
  Description: Performs a database fetch and returns feature objects in
               contig coordinates.
  Returntype : list of Bio::EnsEMBL::*Feature in contig coordinates
  Exceptions : none
  Caller     : BaseFeatureAdaptor, ProxyDnaAlignFeatureAdaptor::generic_fetch

=cut
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
  
sub generic_fetch {
  my ($self, $constraint, $logic_name) = @_;
  
  my $tablename = $self->_tablename();
  my $columns = join(', ', $self->_columns());
  
  if($logic_name) {
    #determine the analysis id via the logic_name
    my $aa = $self->db->get_AnalysisAdaptor();
    my $analysis = $aa->fetch_by_logic_name($logic_name);
    unless(defined $analysis && $analysis->dbID() ) {
      $self->warn("No analysis for logic name $logic_name exists\n");
      return ();
    }
    
    my $analysis_id = $analysis->dbID();
    
    if($constraint) {
      $constraint .= " AND analysis_id = $analysis_id";
    } else {
      $constraint = " analysis_id = $analysis_id";
    }
  } 
      
  my $sql = 
    "SELECT $columns 
     FROM $tablename";

  if($constraint) {
     $sql .= " WHERE $constraint";
  }
  
  my $sth = $self->prepare($sql);
  $sth->execute();
  
  my $hashref;
  my @out;

  while($hashref = $sth->fetchrow_hashref()) {
    push @out, $self->_obj_from_hashref($hashref);
  }
  
  return @out;
}


=head2 fetch_by_dbID
117
118
119
120
121
122
123
124
125

  Arg [1]    : int $id
               the unique database identifier for the feature to be obtained 
  Example    : $feat = $adaptor->fetch_by_dbID(1234);
  Description: Returns the feature created from the database defined by the
               the id $id. 
  Returntype : Bio::EnsEMBL::*Feature
  Exceptions : thrown if $id is not defined
  Caller     : general
126
127
128
129
130
131
132
133
134
135
136
137
138

=cut

sub fetch_by_dbID{
  my ($self,$id) = @_;
  
  unless(defined $id) {
    $self->throw("fetch_by_dbID must have an id");
  }

  my $tablename = $self->_tablename();
  my $constraint = "${tablename}_id = $id";

139
140
141
  #return first element of _generic_fetch list
  my ($feat) = $self->generic_fetch($constraint); 
  return $feat;
142
143
144
}


145
=head2 fetch_by_Contig_constraint
146

147
148
  Arg [1]    : Bio::EnsEMBL::RawContig $contig
               The contig object from which features are to be obtained
149
150
151
152
  Arg [2]    : (optional) string $constraint
               An SQL query constraint (i.e. part of the WHERE clause)
  Arg [3]    : (optional) string $logic_name
               the logic name of the type of features to obtain
153
  Example    : @fts = $a->fetch_by_Contig_constraint($contig, 'perc_ident>5.0');
154
155
156
157
158
159
160
  Description: Returns a list of features created from the database which are 
               are on the contig defined by $cid and fulfill the SQL constraint
               defined by $constraint. If logic name is defined, only features
               with an analysis of type $logic_name will be returned. 
  Returntype : list of Bio::EnsEMBL::*Feature in contig coordinates
  Exceptions : thrown if $cid is not defined
  Caller     : general
161
162
163

=cut

164
165
sub fetch_by_Contig_constraint {
  my ($self, $contig, $constraint, $logic_name) = @_;
166
  
167
168
  unless( defined $contig ) {
    $self->throw("fetch_by_Contig_constraint must have an contig");
169
170
  }

171
172
173
174
175
176
  unless( ref $contig && $contig->isa('Bio::EnsEMBL::RawContig')) {
    $self->throw("contig argument is not a Bio::EnsEMBL::RawContig object\n");
  }

  my $cid = $contig->dbID();

177
178
179
180
181
182
183
184
185
  if($constraint) {
    $constraint .= " AND contig_id = $cid";
  } else {
    $constraint = "contig_id = $cid";
  }

  return $self->generic_fetch($constraint, $logic_name);
}

186

187
=head2 fetch_by_Contig
188

189
190
  Arg [1]    : Bio::EnsEMBL::RawContig $contig 
               the contig from which features should be obtained
191
192
  Arg [2]    : (optional) string $logic_name
               the logic name of the type of features to obtain
193
  Example    : @fts = $a->fetch_by_Contig($contig, 'swall');
194
195
196
197
198
199
200
201
202
  Description: Returns a list of features created from the database which are 
               are on the contig defined by $cid If logic name is defined, 
               only features with an analysis of type $logic_name will be 
               returned. 
  Returntype : list of Bio::EnsEMBL::*Feature in contig coordinates
  Exceptions : none
  Caller     : general

=cut
203
   
204
205
sub fetch_by_Contig{
  my ($self, $contig, $logic_name) = @_;
206
207

  #fetch by contig id constraint with empty constraint
208
  return $self->fetch_by_Contig_constraint($contig, '',$logic_name);
209
210
211
}


212
=head2 fetch_by_Contig_and_score
213

214
215
  Arg [1]    : Bio::EnsEMBL::RawContig $contig 
               the contig from which features should be obtained
216
217
218
219
  Arg [2]    : float $score
               the lower bound of the score of the features to obtain
  Arg [3]    : (optional) string $logic_name
               the logic name of the type of features to obtain
220
  Example    : @fts = $a->fetch_by_Contig_and_score(1, 50.0, 'swall');
221
222
223
224
225
226
227
228
229
230
  Description: Returns a list of features created from the database which are 
               are on the contig defined by $cid and which have score greater
               than score.  If logic name is defined, only features with an 
               analysis of type $logic_name will be returned. 
  Returntype : list of Bio::EnsEMBL::*Feature in contig coordinates
  Exceptions : thrown if $score is not defined
  Caller     : general

=cut

231
232
sub fetch_by_Contig_and_score{
  my($self, $contig, $score, $logic_name) = @_;
233
234
235

  my $constraint;

236
  if(!defined $score){
237
238
239
240
241
242
    $self->throw("need a score even if its 0\n");
  } else{
    $constraint = "score > $score";
  }
    
  my @features = 
243
    $self->fetch_by_Contig_constraint($contig, $constraint, $logic_name);
244
245
246
247
248
  
  return @features;
}


249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
=head2 fetch_by_Slice_constraint

  Arg [1]    : Bio::EnsEMBL::Slice $slice
               the slice from which to obtain features
  Arg [2]    : (optional) string $constraint
               An SQL query constraint (i.e. part of the WHERE clause)
  Arg [3]    : (optional) string $logic_name
               the logic name of the type of features to obtain
  Example    : @fts = $a->fetch_by_Slice_constraint($slice, 'perc_ident > 5');
  Description: Returns a list of features created from the database which are 
               are on the Slice defined by $slice and fulfill the SQL 
               constraint defined by $constraint. If logic name is defined, 
               only features with an analysis of type $logic_name will be 
               returned. 
  Returntype : list of Bio::EnsEMBL::*Feature in Slice coordinates
  Exceptions : thrown if $slice is not defined
  Caller     : Bio::EnsEMBL::Slice

=cut

269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
sub fetch_by_Slice_constraint {
  my($self, $slice, $constraint, $logic_name) = @_;

  if(!$slice){
    $self->throw("need a slice to work\n");
  }
  unless($slice->isa("Bio::EnsEMBL::Slice")) {
    $self->throw("$slice isn't a slice");
  }

  my @features = 
    $self->fetch_by_assembly_location_constraint($slice->chr_start,
						 $slice->chr_end,
						 $slice->chr_name,
						 $slice->assembly_type,
						 $constraint,
						 $logic_name);

  #convert from chromosomal coordinates to slice coordinates
  foreach my $f (@features){
    my $start = ($f->start - ($slice->chr_start - 1));
    my $end = ($f->end - ($slice->chr_start - 1));
   
    $f->start($start);
    $f->end($end);
    $f->attach_seq($slice);
  }

  return @features;
}


301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
=head2 fetch_by_Slice

  Arg [1]    : Bio::EnsEMBL::Slice $slice
               the slice from which to obtain features
  Arg [2]    : (optional) string $logic_name
               the logic name of the type of features to obtain
  Example    : @fts = $a->fetch_by_Slice($slice, 'swall');
  Description: Returns a list of features created from the database which are 
               are on the Slice defined by $slice.If $logic_name is defined, 
               only features with an analysis of type $logic_name will be 
               returned. 
  Returntype : list of Bio::EnsEMBL::*Feature in Slice coordinates
  Exceptions : none
  Caller     : Bio::EnsEMBL::Slice

=cut

318
319
320
321
322
323
324
325
sub fetch_by_Slice {
  my ($self, $slice, $logic_name) = @_;
  
  #fetch by constraint with empty constraint
  return $self->fetch_by_Slice_constraint($slice, '', $logic_name);
}


326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
=head2 fetch_by_Slice_and_score

  Arg [1]    : Bio::EnsEMBL::Slice $slice
               the slice from which to obtain features
  Arg [2]    : float $score
               lower bound of the the score of the features retrieved
  Arg [3]    : (optional) string $logic_name
               the logic name of the type of features to obtain
  Example    : @fts = $a->fetch_by_Slice($slice, 'swall');
  Description: Returns a list of features created from the database which are 
               are on the Slice defined by $slice and which have a score 
               greated than $score. If $logic_name is defined, 
               only features with an analysis of type $logic_name will be 
               returned. 
  Returntype : list of Bio::EnsEMBL::*Feature in Slice coordinates
  Exceptions : none
  Caller     : Bio::EnsEMBL::Slice

=cut

346
347
348
349
sub fetch_by_Slice_and_score {
  my ($self, $slice, $score, $logic_name) = @_;
  my $constraint;

350
  if(!defined $score) {
351
352
353
354
355
356
357
358
359
360
361
    $self->throw("need a score even if its 0\n");
  } else {
    $constraint = "score > $score";
  }

  return $self->fetch_by_Slice_constraint($slice, $constraint, $logic_name);
}  


=head2 fetch_by_assembly_location

362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
  Arg [1]    : int $start
               the start of the assembly region from which to obtain align
               features (in chromosomal coords).
  Arg [2]    : int $end
               the end of the assembly region from which to obtain align 
               features (in chromosomal coords).
  Arg [3]    : string $chr
               the name of the chromosome from which to obtain align features.
  Arg [4]    : string $type
               the type of assembly to obtain features from
  Arg [5]    : (optional) string $logic_name
               the logic name of the type of features to obtain
  Example    : @feats = $adaptor->fetch_by_assembly_location(1, 10000, '9', 'NCBI30');
  Description: Returns a list of features created from the database which are 
               are in the assembly region defined by $start, $end, and $chr. 
               If $logic_name is defined, only features with an analysis 
               of type $logic_name will be returned.
  Returntype : list of Bio::EnsEMBL::*Features
  Exceptions : none
  Caller     : general
382
383
384
385
386
387
388
389
390
391
392
393

=cut

sub fetch_by_assembly_location{
  my ($self,$start,$end,$chr,$type, $logic_name) = @_;

  #fetch by assembly location constraint w/ empty constraint
  return $self->fetch_by_assembly_location_constraint($start, $end, $chr,
						      $type, '', $logic_name);
}


394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
=head2 fetch_by_assembly_location_and_score

  Arg [1]    : int $start
               the start of the assembly region from which to obtain align
               features (in chromosomal coords).
  Arg [2]    : int $end
               the end of the assembly region from which to obtain align 
               features (in chromosomal coords).
  Arg [3]    : string $chr
               the name of the chromosome from which to obtain align features.
  Arg [5]    : string $type
               the type of assembly to obtain features from
  Arg [6]    : float $score
               a lower bound for the score of feats to obtain
  Arg [7]    : (optional) string $logic_name
               the logic name of the type of features to obtain
  Example    : @f = $a->fetch_by_assembly_location_and_score(1,10000,'9','NCBI30', 50.0);
  Description: Returns a list of features created from the database which are 
               are in the assembly region defined by $start, $end, and $chr, 
               and with a percentage identity greater than $pid.  If 
               $logic_name is defined, only features with an analysis of type 
               $logic_name will be returned. 
  Returntype : list of Bio::EnsEMBL::*AlignFeature in chromosomal coordinates
  Exceptions : thrown if $score is not defined
  Caller     : general

=cut

422
423
424
425
sub fetch_by_assembly_location_and_score{
  my ($self, $start, $end, $chr, $type, $score, $logic_name) = @_;
  my $constraint;

426
  if(!defined $score) {
427
428
429
430
431
432
433
434
435
436
437
438
    $self->throw("need a score even if its 0\n");
  } else {
    $constraint = "score > $score";
  }

  return $self->fetch_by_assembly_location_constraint($start, $end, $chr,$type,
						     $constraint, $logic_name);
}


=head2 fetch_by_assembly_location_constraint

439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
  Arg [1]    : int $start
               the start of the assembly region from which to obtain align
               features (in chromosomal coords).
  Arg [2]    : int $end
               the end of the assembly region from which to obtain align 
               features (in chromosomal coords).
  Arg [3]    : string $chr
               the name of the chromosome from which to obtain align features.
  Arg [4]    : string $type
               the type of assembly to obtain features from
  Arg [5]    : (optional) string $constraint
               a SQL constraint restricting which features should be obtained
  Arg [5]    : (optional) string $logic_name
               the logic name of the type of features to obtain
  Example    : @f = $a->fetch_by_assembly_location_constraint(1,10000,'9','NCBI30', 'score > 50.0', 'swall');
  Description: Returns a list of features created from the database which are 
               are in the assembly region defined by $start, $end, and $chr, 
               and with a percentage identity greater than $pid.  If 
               $logic_name is defined, only features with an analysis of type 
               $logic_name will be returned. 
  Returntype : list of Bio::EnsEMBL::*AlignFeature in chromosomal coordinates
  Exceptions : thrown if $score is not defined
  Caller     : BaseFeatureAdaptor
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526

=cut

sub fetch_by_assembly_location_constraint {
  my ($self, $chr_start, $chr_end, $chr, $type, $constraint, $logic_name) = @_;

  if( !defined $type ) {
    $self->throw("Assembly location must be start,end,chr,type");
  }

  if( $chr_start !~ /^\d/ || $chr_end !~ /^\d/ ) {
    $self->throw("start/end must be numbers not $chr_start,$chr_end " .
		 "(have you typed the location in the right way around" .
		 "- start,end,chromosome,type)?");
  }
  
  my $mapper = $self->db->get_AssemblyMapperAdaptor->fetch_by_type($type);  
  my @cids = $mapper->list_contig_ids($chr, $chr_start ,$chr_end);
  
  if( scalar(@cids) == 0 ) {
    return ();
  }

  my $cid_list = join(',',@cids);

  #construct the constraint for the contig ids 
  if($constraint) {
    $constraint .= " AND contig_id IN ($cid_list)";
  } else {
    $constraint = "contig_id IN ($cid_list)";
  }
  
  my @features = $self->generic_fetch($constraint, $logic_name); 

  my @out;

  #convert the features to assembly coordinates from raw contig coordinates
  foreach my $f (@features) {
    #since feats were obtained in contig coords, attached seq is a contig
    my $contig_id = $f->entire_seq->dbID();
    my @coord_list = 
      $mapper->map_coordinates_to_assembly($contig_id, $f->start(),
					   $f->end(),$f->strand(),"rawcontig");
       
    # coord list > 1 - means does not cleanly map At the moment, skip
    if( scalar(@coord_list) > 1 ) {
      next;
    }
     
    my ($coord) = @coord_list;

    #maps to a gap in assembly?
    if($coord->isa("Bio::EnsEMBL::Mapper::Gap")){
      next;
    }

    #maps to region outside desired area
    if(!($coord->start() >= $chr_start) ||
       !($coord->end() <= $chr_end)) {
      next;
    }
        
    $f->start($coord->start());
    $f->end($coord->end());
    $f->strand($coord->strand());
527
528
529
530
531
    #$f->seqname($coord->id());

    #
    # Should we attach a slice of the entire chromosome here? (mcvicker)
    #
532
533
534
535
536
537
538
539

    push(@out,$f);
  }
  
  return @out;
}


540
=head2 store
541

542
543
544
545
546
547
548
549
  Arg [1]    : list of Bio::EnsEMBL::*Feature
  Example    : $adaptor->store(@feats);
  Description: ABSTRACT  Subclasses are responsible for implementing this 
               method.  It should take a list of features and store them in 
               the database.
  Returntype : none
  Exceptions : thrown method is not implemented by subclass
  Caller     : general
550
551
552
553
554
555
556
557
558
559
560
561

=cut

sub store{
  my $self = @_;

  $self->throw("Abstract method store not defined by implementing subclass\n");
}


=head2 _tablename

562
563
564
565
566
567
568
569
  Args       : none
  Example    : $tablename = $self->_table_name()
  Description: ABSTRACT PROTECTED Subclasses are responsible for implementing
               this method.  It should return the name of the table to be
               used to obtain features.  
  Returntype : string
  Exceptions : thrown if not implemented by subclass
  Caller     : BaseFeatureAdaptor::generic_fetch
570
571
572
573
574
575
576
577
578
579
580

=cut

sub _tablename {
  my $self = shift;

  $self->throw("abstract method _tablename not defined by implementing" .
               " subclass of AlignFeatureAdaptor");
  return undef;
}

581

582
583
=head2 _columns

584
585
586
587
588
589
590
591
  Args       : none
  Example    : $tablename = $self->_columns()
  Description: ABSTRACT PROTECTED Subclasses are responsible for implementing
               this method.  It should return a list of columns to be used
               for feature creation
  Returntype : list of strings
  Exceptions : thrown if not implemented by subclass
  Caller     : BaseFeatureAdaptor::generic_fetch
592
593
594
595
596
597
598
599
600
601
602

=cut

sub _columns {
  my $self = shift;

  $self->throw("abstract method _columns not defined by implementing" .
               " subclass of AlignFeatureAdaptor");
}


603
604
605
606
607
608
609
610
611
612
613
614
615
616
=head2 _obj_from_hashref

  Arg [1]    : DBI::row_hashref $hashref containing key-value pairs 
               for each of the columns specified by the _columns method
  Example    : my @feats = $self->_obj_from_hashref
  Description: ABSTRACT PROTECTED The subclass is responsible for implementing
               this method.  It should take in a DBI row hash reference and
               return a list of created features in contig coordinates.
  Returntype : list of Bio::EnsEMBL::*Features in contig coordinates
  Exceptions : thrown if not implemented by subclass
  Caller     : BaseFeatureAdaptor::generic_fetch

=cut

617
618
619
620
621
622
623
624
625
626
sub _obj_from_hashref {
  my $self = shift;

  $self->throw("abstract method _obj_from_hashref not defined by implementing"
             . " subclass of AlignFeatureAdaptor");
} 

1;