branch-1-6/Bio/Tools/Grail.pm

   1 # $Id$
   2 #
   3 # BioPerl module for Bio::Tools::Grail
   4 #
   5 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   6 #
   7 # Cared for by Jason Stajich <jason-at-bioperl.org>
   8 #
   9 # Copyright Jason Stajich
  10 #
  11 # You may distribute this module under the same terms as perl itself
  12
  13 # POD documentation - main docs before the code
  14
  15 =head1 NAME
  16
  17 Bio::Tools::Grail - Results of one Grail run
  18
  19 =head1 SYNOPSIS
  20
  21    $grail = Bio::Tools::Grail->new(-file => 'result.grail');
  22    # filehandle:
  23    $grail = Bio::Tools::Grail->new( -fh  => \*INPUT );
  24
  25    # parse the results
  26    while($gene = $grail->next_prediction()) {
  27        # $gene is an instance of Bio::Tools::Prediction::Gene
  28
  29        # $gene->exons() returns an array of
  30        # Bio::Tools::Prediction::Exon objects
  31        # all exons:
  32        @exon_arr = $gene->exons();
  33
  34        # initial exons only
  35        @init_exons = $gene->exons('Initial');
  36        # internal exons only
  37        @intrl_exons = $gene->exons('Internal');
  38        # terminal exons only
  39        @term_exons = $gene->exons('Terminal');
  40        # singleton exons only -- should be same as $gene->exons() because
  41        # there are no other exons supposed to exist in this structure
  42        @single_exons = $gene->exons('Single');
  43    }
  44
  45    # essential if you gave a filename at initialization (otherwise the file
  46    # will stay open)
  47    $genscan->close();
  48
  49 =head1 DESCRIPTION
  50
  51 The Grail module provides a parser for Grail gene structure prediction
  52 output.
  53
  54
  55 =head1 FEEDBACK
  56
  57 =head2 Mailing Lists
  58
  59 User feedback is an integral part of the evolution of this and other
  60 Bioperl modules. Send your comments and suggestions preferably to one
  61 of the Bioperl mailing lists.  Your participation is much appreciated.
  62
  63   bioperl-l@bioperl.org                  - General discussion
  64   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  65
  66 =head2 Support
  67
  68 Please direct usage questions or support issues to the mailing list:
  69
  70 I<bioperl-l@bioperl.org>
  71
  72 rather than to the module maintainer directly. Many experienced and
  73 reponsive experts will be able look at the problem and quickly
  74 address it. Please include a thorough description of the problem
  75 with code and data examples if at all possible.
  76
  77 =head2 Reporting Bugs
  78
  79 Report bugs to the Bioperl bug tracking system to help us keep track
  80 the bugs and their resolution.  Bug reports can be submitted via the
  81 web:
  82
  83   http://bugzilla.open-bio.org/
  84
  85 =head1 AUTHOR - Jason Stajich
  86
  87 Email jason-at-bioperl.org
  88
  89 =head1 APPENDIX
  90
  91 The rest of the documentation details each of the object
  92 methods. Internal methods are usually preceded with a _
  93
  94 =cut
  95
  96 # Let the code begin...
  97
  98 package Bio::Tools::Grail;
  99 use strict;
 100
 101 use Bio::Tools::Prediction::Gene;
 102 use Bio::Tools::Prediction::Exon;
 103 use Symbol;
 104
 105 use base qw(Bio::Root::IO Bio::Root::Root);
 106
 107 sub new {
 108   my($class,@args) = @_;
 109
 110   my $self = $class->SUPER::new(@args);
 111   $self->_initialize_io(@args);
 112
 113   return $self;
 114 }
 115
 116 =head2 next_prediction
 117
 118  Title   : next_prediction
 119  Usage   : while($gene = $grail->next_prediction()) {
 120                   # do something
 121            }
 122  Function: Returns the next gene structure prediction of the Grail result
 123            file. Call this method repeatedly until FALSE is returned.
 124
 125  Example :
 126  Returns : A Bio::Tools::Prediction::Gene object.
 127  Args    :
 128
 129 =cut
 130
 131 sub next_prediction {
 132     my ($self) = @_;
 133
 134     # get next gene structure
 135     my $gene = $self->_prediction();
 136
 137     if($gene) {
 138         # fill in predicted protein, and if available the predicted CDS
 139         #
 140         my ($id, $seq);
 141         # use the seq stack if there's a seq on it
 142         my $seqobj = pop(@{$self->{'_seqstack'}});
 143         if(! $seqobj) {
 144             # otherwise read from input stream
 145             ($id, $seq) = $self->_read_fasta_seq();
 146             $seqobj = Bio::PrimarySeq->new('-seq' => $seq,
 147                                            '-display_id' => $id,
 148                                            '-alphabet' => "protein");
 149         }
 150         # check that prediction number matches the prediction number
 151         # indicated in the sequence id (there may be incomplete gene
 152         # predictions that contain only signals with no associated protein
 153         # and CDS, like promoters, poly-A sites etc)
 154         $gene->primary_tag() =~ /[^0-9]([0-9]+)$/;
 155         my $prednr = $1;
 156         if($seqobj->display_id() !~ /_predicted_\w+_$prednr\|/) {
 157             # this is not our sequence, so push back for the next prediction
 158             push(@{$self->{'_seqstack'}}, $seqobj);
 159         } else {
 160             $gene->predicted_protein($seqobj);
 161             # CDS prediction, too?
 162             if($self->_has_cds()) {
 163                 ($id, $seq) = $self->_read_fasta_seq();
 164                 $seqobj = Bio::PrimarySeq->new('-seq' => $seq,
 165                                                '-display_id' => $id,
 166                                                '-alphabet' => "dna");
 167                 $gene->predicted_cds($seqobj);
 168             }
 169         }
 170     }
 171     return $gene;
 172 }
 173
 174 =head2 _parse_predictions
 175
 176  Title   : _parse_predictions()
 177  Usage   : $obj->_parse_predictions()
 178  Function: Parses the prediction section. Automatically called by
 179            next_prediction() if not yet done.
 180  Example :
 181  Returns :
 182
 183 =cut
 184
 185 sub _parse_predictions {
 186     my ($self) = @_;
 187
 188     # code needs to go here
 189
 190     $self->_predictions_parsed(1);
 191 }
 192
 193 =head2 _prediction
 194
 195  Title   : _prediction()
 196  Usage   : $gene = $obj->_prediction()
 197  Function: internal
 198  Example :
 199  Returns :
 200
 201 =cut
 202
 203 sub _prediction {
 204     my ($self) = @_;
 205
 206     return unless(exists($self->{'_preds'}) && @{$self->{'_preds'}});
 207     return shift(@{$self->{'_preds'}});
 208 }
 209
 210 =head2 _add_prediction
 211
 212  Title   : _add_prediction()
 213  Usage   : $obj->_add_prediction($gene)
 214  Function: internal
 215  Example :
 216  Returns :
 217
 218 =cut
 219
 220 sub _add_prediction {
 221     my ($self, $gene) = @_;
 222
 223     if(! exists($self->{'_preds'})) {
 224         $self->{'_preds'} = [];
 225     }
 226     push(@{$self->{'_preds'}}, $gene);
 227 }
 228
 229 =head2 _predictions_parsed
 230
 231  Title   : _predictions_parsed
 232  Usage   : $obj->_predictions_parsed
 233  Function: internal
 234  Example :
 235  Returns : TRUE or FALSE
 236
 237 =cut
 238
 239 sub _predictions_parsed {
 240     my ($self, $val) = @_;
 241
 242     $self->{'_preds_parsed'} = $val if $val;
 243     if(! exists($self->{'_preds_parsed'})) {
 244         $self->{'_preds_parsed'} = 0;
 245     }
 246     return $self->{'_preds_parsed'};
 247 }
 248
 249 =head2 _has_cds
 250
 251  Title   : _has_cds()
 252  Usage   : $obj->_has_cds()
 253  Function: Whether or not the result contains the predicted CDSs, too.
 254  Example :
 255  Returns : TRUE or FALSE
 256
 257 =cut
 258
 259 sub _has_cds {
 260     my ($self, $val) = @_;
 261
 262     $self->{'_has_cds'} = $val if $val;
 263     if(! exists($self->{'_has_cds'})) {
 264         $self->{'_has_cds'} = 0;
 265     }
 266     return $self->{'_has_cds'};
 267 }
 268
 269 1;