|
0
|
1 # $Id: PairwiseStatistics.pm,v 1.2 2002/10/22 07:45:10 lapp Exp $
|
|
|
2 #
|
|
|
3 # BioPerl module for Bio::Align::PairwiseStatistics
|
|
|
4 #
|
|
|
5 # Cared for by Jason Stajich <jason@bioperl.org>
|
|
|
6 #
|
|
|
7 # Copyright Jason Stajich
|
|
|
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::Align::PairwiseStatistics - Base statistic object for Pairwise Alignments
|
|
|
16
|
|
|
17 =head1 SYNOPSIS
|
|
|
18
|
|
|
19 Give standard usage here
|
|
|
20
|
|
|
21 =head1 DESCRIPTION
|
|
|
22
|
|
|
23 Describe the object here
|
|
|
24
|
|
|
25 =head1 FEEDBACK
|
|
|
26
|
|
|
27 =head2 Mailing Lists
|
|
|
28
|
|
|
29 User feedback is an integral part of the evolution of this and other
|
|
|
30 Bioperl modules. Send your comments and suggestions preferably to
|
|
|
31 the Bioperl mailing list. Your participation is much appreciated.
|
|
|
32
|
|
|
33 bioperl-l@bioperl.org - General discussion
|
|
|
34 http://bioperl.org/MailList.shtml - About the mailing lists
|
|
|
35
|
|
|
36 =head2 Reporting Bugs
|
|
|
37
|
|
|
38 Report bugs to the Bioperl bug tracking system to help us keep track
|
|
|
39 of the bugs and their resolution. Bug reports can be submitted via
|
|
|
40 email or the web:
|
|
|
41
|
|
|
42 bioperl-bugs@bioperl.org
|
|
|
43 http://bugzilla.bioperl.org/
|
|
|
44
|
|
|
45 =head1 AUTHOR - Jason Stajich
|
|
|
46
|
|
|
47 Email jason@bioperl.org
|
|
|
48
|
|
|
49 Describe contact details here
|
|
|
50
|
|
|
51 =head1 CONTRIBUTORS
|
|
|
52
|
|
|
53 Additional contributors names and emails here
|
|
|
54
|
|
|
55 =head1 APPENDIX
|
|
|
56
|
|
|
57 The rest of the documentation details each of the object methods.
|
|
|
58 Internal methods are usually preceded with a _
|
|
|
59
|
|
|
60 =cut
|
|
|
61
|
|
|
62
|
|
|
63 # Let the code begin...
|
|
|
64
|
|
|
65
|
|
|
66 package Bio::Align::PairwiseStatistics;
|
|
|
67 use vars qw(@ISA $GapChars);
|
|
|
68 use strict;
|
|
|
69
|
|
|
70 use Bio::Align::StatisticsI;
|
|
|
71 use Bio::Root::Root;
|
|
|
72
|
|
|
73 BEGIN { $GapChars = '(\.|\-)'; }
|
|
|
74
|
|
|
75 @ISA = qw(Bio::Root::Root Bio::Align::StatisticsI );
|
|
|
76
|
|
|
77 =head2 number_of_comparable_bases
|
|
|
78
|
|
|
79 Title : number_of_comparable_bases
|
|
|
80 Usage : my $bases = $stat->number_of_comparable_bases($aln);
|
|
|
81 Function: Returns the count of the number of bases that can be
|
|
|
82 compared (L) in this alignment ( length - gaps)
|
|
|
83 Returns : integer
|
|
|
84 Args : Bio::Align::AlignI
|
|
|
85
|
|
|
86
|
|
|
87 =cut
|
|
|
88
|
|
|
89 sub number_of_comparable_bases{
|
|
|
90 my ($self,$aln) = @_;
|
|
|
91 if( ! defined $aln || ! $aln->isa('Bio::Align::AlignI') ) {
|
|
|
92 $self->warn("Must provide a Bio::Align::AlignI compliant object to Bio::Align::PairwiseStatistics");
|
|
|
93 return 0;
|
|
|
94 } elsif( $aln->no_sequences != 2 ) {
|
|
|
95 $self->warn("only pairwise calculations currently supported");
|
|
|
96 }
|
|
|
97 my $L = $aln->length - $self->number_of_gaps($aln);
|
|
|
98 return $L;
|
|
|
99 }
|
|
|
100
|
|
|
101 =head2 number_of_differences
|
|
|
102
|
|
|
103 Title : number_of_differences
|
|
|
104 Usage : my $nd = $stat->number_of_distances($aln);
|
|
|
105 Function: Returns the number of differences between two
|
|
|
106 Returns : integer
|
|
|
107 Args : Bio::Align::AlignI
|
|
|
108
|
|
|
109
|
|
|
110 =cut
|
|
|
111
|
|
|
112 sub number_of_differences{
|
|
|
113 my ($self,$aln) = @_;
|
|
|
114 if( ! defined $aln || ! $aln->isa('Bio::Align::AlignI') ) {
|
|
|
115 $self->warn("Must provide a Bio::Align::AlignI compliant object to Bio::Align::PairwiseStatistics");
|
|
|
116 return 0;
|
|
|
117 } elsif( $aln->no_sequences != 2 ) {
|
|
|
118 $self->warn("only pairwise calculations currently supported");
|
|
|
119 }
|
|
|
120 my (@seqs);
|
|
|
121 foreach my $seq ( $aln->each_seq) {
|
|
|
122 push @seqs, [ split(//,$seq->seq())];
|
|
|
123 }
|
|
|
124 my $firstseq = shift @seqs;
|
|
|
125 # my $secondseq = shift @seqs;
|
|
|
126 my $diffcount = 0;
|
|
|
127 for (my $i = 0;$i<$aln->length; $i++ ) {
|
|
|
128 next if( $firstseq->[$i] =~ /^$GapChars$/);
|
|
|
129 foreach my $seq ( @seqs ) {
|
|
|
130 next if( $seq->[$i] =~ /^$GapChars$/);
|
|
|
131 if( $firstseq->[$i] ne $seq->[$i] ) {
|
|
|
132 $diffcount++;
|
|
|
133 }
|
|
|
134 }
|
|
|
135 }
|
|
|
136 return $diffcount;
|
|
|
137 }
|
|
|
138
|
|
|
139 =head2 number_of_gaps
|
|
|
140
|
|
|
141 Title : number_of_gaps
|
|
|
142 Usage : my $nd = $stat->number_of_gaps($aln);
|
|
|
143 Function: Returns the number of differences between two
|
|
|
144 Example :
|
|
|
145 Returns :
|
|
|
146 Args :
|
|
|
147
|
|
|
148
|
|
|
149 =cut
|
|
|
150
|
|
|
151 sub number_of_gaps{
|
|
|
152 my ($self,$aln) = @_;
|
|
|
153 if( ! defined $aln || ! $aln->isa('Bio::Align::AlignI') ) {
|
|
|
154 $self->warn("Must provide a Bio::Align::AlignI compliant object to Bio::Align::PairwiseStatistics");
|
|
|
155 return 0;
|
|
|
156 } elsif( $aln->no_sequences != 2 ) {
|
|
|
157 $self->warn("only pairwise calculations currently supported");
|
|
|
158 }
|
|
|
159 my (@seqs);
|
|
|
160 foreach my $seq ( $aln->each_seq) {
|
|
|
161 push @seqs, [ split(//,$seq->seq())];
|
|
|
162 }
|
|
|
163 my $firstseq = shift @seqs;
|
|
|
164 # my $secondseq = shift @seqs;
|
|
|
165 my $gapcount = 0;
|
|
|
166 for (my $i = 0;$i<$aln->length; $i++ ) {
|
|
|
167 ($gapcount++) && next if( $firstseq->[$i] =~ /^$GapChars$/);
|
|
|
168 foreach my $seq ( @seqs ) {
|
|
|
169 ($gapcount++) && next if( $seq->[$i] =~ /^$GapChars$/);
|
|
|
170 }
|
|
|
171 }
|
|
|
172 return $gapcount;
|
|
|
173 }
|
|
|
174
|
|
|
175 1;
|
