File Coverage

blib/lib/Bio/MUST/Core/SeqMask/Freqs.pm
Criterion Covered Total %
statement 61 61 100.0
branch 2 2 100.0
condition 3 4 75.0
subroutine 10 10 100.0
pod 1 1 100.0
total 77 78 98.7


line stmt bran cond sub pod time code
1             package Bio::MUST::Core::SeqMask::Freqs;
2             # ABSTRACT: Arbitrary frequencies for sequence sites
3             $Bio::MUST::Core::SeqMask::Freqs::VERSION = '0.212670';
4 17     17   11170 use Moose;
  17         54  
  17         116  
5 17     17   120671 use namespace::autoclean;
  17         50  
  17         148  
6              
7 17     17   1381 use autodie;
  17         45  
  17         121  
8 17     17   91633 use feature qw(say);
  17         50  
  17         1374  
9              
10 17     17   127 use List::AllUtils qw(sum);
  17         54  
  17         1152  
11              
12 17     17   133 use Bio::MUST::Core::Types;
  17         40  
  17         579  
13 17     17   126 use aliased 'Bio::MUST::Core::SeqMask::Rates';
  17         41  
  17         135  
14              
15              
16             # public hash containing freqs by sequence and site
17             # Note: this hash is actually a Tie::IxHash (see Profiles)
18             has 'freq_for_at' => (
19             traits => ['Hash'],
20             is => 'ro',
21             isa => 'HashRef[ArrayRef[Num]]',
22             required => 1,
23             handles => {
24             count_freqs_at => 'count',
25             all_freqs_at => 'values',
26             freqs_at_for => 'get',
27             all_ids => 'keys',
28             },
29             );
30              
31              
32             # private SeqMask::Rates-like object derived by averaging freqs over seqs
33             has '_mask' => (
34             is => 'ro',
35             isa => 'Bio::MUST::Core::SeqMask::Rates',
36             init_arg => undef,
37             lazy => 1,
38             builder => '_build_mask',
39             handles => {
40             mask_len => 'mask_len',
41             all_freqs => 'all_states',
42             min_freq => 'min_rate',
43             max_freq => 'max_rate',
44             bin_freqs_masks => 'bin_rates_masks',
45             freqs_mask => 'rates_mask',
46             },
47             );
48              
49              
50             # private hash containing freqs averaged over sites
51             has '_avg_freq_for' => (
52             traits => ['Hash'],
53             is => 'ro',
54             isa => 'HashRef[Num]',
55             init_arg => undef,
56             lazy => 1,
57             builder => '_build_avg_freq_for',
58             handles => {
59             avg_freq_for => 'get',
60             },
61             );
62              
63              
64             ## no critic (ProhibitUnusedPrivateSubroutines)
65              
66             sub _build_mask {
67 3     3   7 my $self = shift;
68              
69 3         7 my @mask;
70              
71             # average freqs over seqs
72 3         125 for my $freqs_at ($self->all_freqs_at) {
73 161         1358 my $i = 0;
74 161         208 $mask[$i++] += $_ for @{$freqs_at};
  161         5470  
75             }
76 3         216 my $n = $self->count_freqs_at;
77 3         520 @mask = map { $_ / $n } @mask;
  555         872  
78              
79 3         144 return Rates->new( mask => \@mask );
80             }
81              
82              
83             sub _build_avg_freq_for {
84 2     2   5 my $self = shift;
85              
86 2         5 my %avg_freq_for;
87              
88             # average freqs over sites
89 2         13 my $n = $self->mask_len;
90 2         78 for my $id ($self->all_ids) {
91 160         2310 $avg_freq_for{$id} = sum( @{ $self->freqs_at_for($id) } ) / $n;
  160         6044  
92             }
93              
94 2         122 return \%avg_freq_for;
95             }
96              
97             ## use critic
98              
99              
100              
101             sub store {
102 3     3 1 12 my $self = shift;
103 3         9 my $outfile = shift;
104 3   50     19 my $args = shift // {}; # HashRef (should not be empty...)
105              
106 3   100     47 my $reorder = $args->{reorder} // 0;
107              
108 3         21 open my $out, '>', $outfile;
109              
110             # optionally sort ids by descending average freq (over sites)
111 3         3005 my @ids = $self->all_ids;
112             @ids = sort {
113 3 100       977 $self->avg_freq_for($b) <=> $self->avg_freq_for($a)
  411         15171  
114             } @ids if $reorder;
115              
116             # output header
117 3         7 say {$out} join "\t", 'site', 'f(i,.)', @ids;
  3         51  
118              
119             # output average freqs (over sites)
120 3         8 say {$out} join "\t", 'f(.,j)', q{}, map { $self->avg_freq_for($_) } @ids;
  3         23  
  240         9016  
121              
122             # setup rows with site numbers and average freqs (over seqs)
123 3         33 my @rows = 1..$self->mask_len;
124 3         18 my @avg_freqs_at = $self->all_freqs;
125 3         937 $_ .= "\t" . shift @avg_freqs_at for @rows;
126              
127             # assemble freqs by site (one site by row)
128 3         10 for my $id (@ids) {
129 240         441 my @freqs_at = @{ $self->freqs_at_for($id) };
  240         9951  
130 240         55555 $_ .= "\t" . shift @freqs_at for @rows;
131             }
132              
133             # output freqs for all sites
134 3         27 say {$out} $_ for @rows;
  555         2164  
135              
136 3         222 return;
137             }
138              
139             __PACKAGE__->meta->make_immutable;
140             1;
141              
142             __END__
143              
144             =pod
145              
146             =head1 NAME
147              
148             Bio::MUST::Core::SeqMask::Freqs - Arbitrary frequencies for sequence sites
149              
150             =head1 VERSION
151              
152             version 0.212670
153              
154             =head1 SYNOPSIS
155              
156             # TODO
157              
158             =head1 DESCRIPTION
159              
160             # TODO
161              
162             =head1 METHODS
163              
164             =head2 store
165              
166             =head1 AUTHOR
167              
168             Denis BAURAIN <denis.baurain@uliege.be>
169              
170             =head1 COPYRIGHT AND LICENSE
171              
172             This software is copyright (c) 2013 by University of Liege / Unit of Eukaryotic Phylogenomics / Denis BAURAIN.
173              
174             This is free software; you can redistribute it and/or modify it under
175             the same terms as the Perl 5 programming language system itself.
176              
177             =cut