File Coverage

Bio/Tools/Phylo/Gerp.pm
Criterion Covered Total %
statement 27 27 100.0
branch 1 2 50.0
condition 3 4 75.0
subroutine 6 6 100.0
pod 2 2 100.0
total 39 41 95.1


line stmt bran cond sub pod time code
1             # $Id: Gumby.pm,v 1.2 2007/06/14 18:01:52 nathan Exp $
2             #
3             # BioPerl module for Bio::Tools::Phylo::Gerp
4             #
5             # Please direct questions and support issues to
6             #
7             # Cared for by Sendu Bala
8             #
9             # Copyright Sendu Bala
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::Phylo::Gerp - Parses output from GERP
18              
19             =head1 SYNOPSIS
20              
21             use strict;
22              
23             use Bio::Tools::Phylo::Gerp;
24              
25             my $parser = Bio::Tools::Phylo::Gerp->new(-file => "alignment.rates.elems");
26              
27             while (my $feat = $parser->next_result) {
28             my $start = $feat->start;
29             my $end = $feat->end;
30             my $rs_score = $feat->score;
31             my $p_value = ($feat->annotation->get_Annotations('p-value'))[0]->value;
32             }
33              
34             =head1 DESCRIPTION
35              
36             This module is used to parse the output from 'GERP' (v2) by Eugene Davydov
37             (originally Gregory M. Cooper et al.). You can get details here:
38             http://mendel.stanford.edu/sidowlab/
39              
40             It works on the .elems files produced by gerpelem.
41              
42             Each result is a Bio::SeqFeature::Annotated representing a single constrained
43             element.
44              
45             =head1 FEEDBACK
46              
47             =head2 Mailing Lists
48              
49             User feedback is an integral part of the evolution of this and other
50             Bioperl modules. Send your comments and suggestions preferably to
51             the Bioperl mailing list. Your participation is much appreciated.
52              
53             bioperl-l@bioperl.org - General discussion
54             http://bioperl.org/wiki/Mailing_lists - About the mailing lists
55              
56             =head2 Support
57              
58             Please direct usage questions or support issues to the mailing list:
59              
60             I
61              
62             rather than to the module maintainer directly. Many experienced and
63             reponsive experts will be able look at the problem and quickly
64             address it. Please include a thorough description of the problem
65             with code and data examples if at all possible.
66              
67             =head2 Reporting Bugs
68              
69             Report bugs to the Bioperl bug tracking system to help us keep track
70             of the bugs and their resolution. Bug reports can be submitted via the
71             web:
72              
73             https://github.com/bioperl/bioperl-live/issues
74              
75             =head1 AUTHOR - Sendu Bala
76              
77             Email bix@sendu.me.uk
78              
79             =head1 APPENDIX
80              
81             The rest of the documentation details each of the object methods.
82             Internal methods are usually preceded with a _
83              
84             =cut
85              
86             # Let the code begin...
87              
88             package Bio::Tools::Phylo::Gerp;
89 1     1   479 use strict;
  1         1  
  1         25  
90              
91 1     1   326 use Bio::SeqFeature::Generic;
  1         2  
  1         25  
92 1     1   6 use Bio::Annotation::SimpleValue;
  1         1  
  1         20  
93              
94 1     1   3 use base qw(Bio::Root::Root Bio::Root::IO);
  1         1  
  1         271  
95              
96              
97             =head2 new
98              
99             Title : new
100             Usage : my $obj = Bio::Tools::Phylo::Gerp->new();
101             Function: Builds a new Bio::Tools::Phylo::Gerp object
102             Returns : Bio::Tools::Phylo::Gerp
103             Args : -file (or -fh) should contain the contents of a gerpelem .elems file
104              
105             =cut
106              
107             sub new {
108 1     1 1 2 my ($class, @args) = @_;
109 1         9 my $self = $class->SUPER::new(@args);
110            
111 1         6 $self->_initialize_io(@args);
112            
113 1         7 return $self;
114             }
115              
116             =head2 next_result
117              
118             Title : next_result
119             Usage : $result = $obj->next_result();
120             Function: Returns the next result available from the input, or undef if there
121             are no more results.
122             Returns : Bio::SeqFeature::Annotated object. Features are annotated with a tag
123             for 'pvalue', and a 'predicted' tag. They have no sequence id unless
124             the input GERP file is non-standard, with the seq id as the 6th
125             column.
126              
127             NB: feature coordinates are alignment columns of the alignment
128             used to create the result file.
129             Args : none
130              
131             =cut
132              
133             sub next_result {
134 6     6 1 6 my ($self) = @_;
135            
136 6   100     17 my $line = $self->_readline || return;
137            
138 5         29 while ($line !~ /^\d+\s+\d+\s+\d+\s+\S+\s+\S+\s*(?:\S+\s*)?$/) {
139 5   50     7 $line = $self->_readline || return;
140             }
141            
142             #start end length RS-score p-value
143             # code elsewhere adds seq_id on the end (not valid GERP), so we capture that
144             # if present
145 5         26 my ($start, $end, undef, $rs_score, $p_value, $seq_id) = split(/\s+/, $line);
146 5 50       29 my $feat = Bio::SeqFeature::Generic->new(
147             $seq_id ? (-seq_id => $seq_id) : (),
148             -start => $start,
149             -end => $end,
150             -strand => 1,
151             -score => $rs_score,
152             #-type => 'conserved_region', ***causes 740x increase in SeqFeatureDB storage requirments!
153             -source => 'GERP');
154            
155 5         20 my $sv = Bio::Annotation::SimpleValue->new(-tagname => 'predicted', -value => 1);
156 5         12 $feat->annotation->add_Annotation($sv);
157 5         12 $sv = Bio::Annotation::SimpleValue->new(-tagname => 'pvalue', -value => $p_value);
158 5         10 $feat->annotation->add_Annotation($sv);
159            
160 5         10 return $feat;
161             }
162              
163             1;