File Coverage

lib/Statistics/Sampler/Multinomial.pm
Criterion Covered Total %
statement 72 81 88.8
branch 12 18 66.6
condition 3 9 33.3
subroutine 15 18 83.3
pod 4 4 100.0
total 106 130 81.5


line stmt bran cond sub pod time code
1             package Statistics::Sampler::Multinomial;
2            
3 2     2   89324 use 5.014;
  2         4  
4 2     2   7 use warnings;
  2         2  
  2         42  
5 2     2   6 use strict;
  2         2  
  2         50  
6            
7             our $VERSION = '0.5';
8            
9 2     2   6 use Carp;
  2         3  
  2         98  
10 2     2   440 use Ref::Util qw /is_arrayref/;
  2         431  
  2         100  
11 2     2   23 use List::Util qw /min sum/;
  2         3  
  2         129  
12 2     2   503 use List::MoreUtils qw /first_index/;
  2         6581  
  2         8  
13 2     2   693 use Scalar::Util qw /blessed/;
  2         4  
  2         929  
14            
15             sub new {
16 7     7 1 38510 my ($class, %args) = @_;
17            
18 7         9 my $data = $args{data};
19 7 100       31 croak 'data argument not passed'
20             if !defined $data;
21 6 100       32 croak 'data argument is not an array ref'
22             if !is_arrayref ($data);
23            
24 3     35   27 my $first_neg_idx = first_index {$_ < 0} @$data;
  35         24  
25 3 100       21 croak "negative values passed in data array"
26             if $first_neg_idx >= 0;
27            
28            
29             my $self = {
30             data => $data,
31             data_sum_to_one => $args{data_sum_to_one},
32 2         6 };
33            
34 2         3 bless $self, $class;
35            
36 2         2 my $prng = $args{prng};
37 2         5 $self->_validate_prng_object ($prng);
38             $self->{prng}
39 2   33     56 = $prng
40             // "${class}::DefaultPRNG"->new;
41            
42 2         23 return $self;
43             }
44            
45             sub _validate_prng_object {
46 2     2   2 my ($self, $prng) = @_;
47            
48             # Math::Random::MT::Auto has boolean op overloading
49             # so make sure we don't trigger it or our tests fail
50             # i.e. don't use "if $prng" directly
51             # (and we waste a random number, but that's less of an issue)
52             #return 1 if !defined $prng;
53            
54             # no default yet, so croak if not passed
55 2 50       3 croak "prng arg is not defined. "
56             . "Math::Random::MT::Auto is a good default to use"
57             if not defined $prng;
58 2 50       8 croak 'prng arg is not an object'
59             if not blessed $prng;
60 2 50       10 croak 'prng arg does not have binomial() method'
61             if not $prng->can('binomial');
62            
63 2         21 return 1;
64             }
65            
66            
67             sub _initialise {
68 2     2   3 my ($self, %args) = @_;
69            
70 2         15 my $probs = $self->{data};
71            
72             # caller has not promised they sum to 1
73 2 50       3 if ( $self->{data_sum_to_one}) {
74 0         0 $self->{sum} = 1;
75             }
76             else {
77 2         11 my $sum = sum (@$probs);
78 2 50       5 if ($sum != 1) {
79 2         3 my @scaled_probs = map {$_ / $sum} @$probs;
  34         28  
80 2         3 $probs = \@scaled_probs;
81             }
82 2         4 $self->{sum} = $sum;
83             }
84            
85 2         4 return;
86             }
87            
88             sub get_class_count {
89 0     0 1 0 my $self = shift;
90 0         0 my $aref = $self->{data};
91 0         0 return scalar @$aref;
92             }
93            
94             sub draw {
95 1     1 1 161 my ($self, $args) = @_;
96            
97             # inefficient, but why would one want a single draw?
98             # here for compatibility with SSM::AliasMethod
99 1         4 return $self->draw_n_samples (1);
100             }
101            
102             sub draw_n_samples {
103 2     2 1 7 my ($self, $n) = @_;
104            
105 2         3 my $prng = $self->{prng};
106            
107             my $data = $self->{data}
108 2   33     5 // croak 'it appears setup has not been run yet';
109 2         3 my $K = scalar @$data;
110 2   33     6 my $norm = $self->{sum} // do {$self->_initialise; $self->{sum}};
  2         3  
  2         2  
111            
112 2         1 my @draws;
113 2         14 my ($sum_p, $sum_n) = (0, 0);
114            
115 2         6 foreach my $kk (0..($K-1)) {
116 34 50       38 if (!$data->[$kk]) {
117 0         0 $draws[$kk] = 0;
118 0         0 next;
119             }
120            
121 34         76 my $res = $prng->binomial (
122             # MRMA does not like p>1
123             min (1, $data->[$kk] / ($norm - $sum_p)),
124             ($n - $sum_n),
125             );
126 34         27 $draws[$kk] = $res;
127 34         20 $sum_p += $data->[$kk];
128 34         24 $sum_n += $res;
129             }
130            
131 2         5 return \@draws;
132             }
133            
134             # Cuckoo package to act as a method wrapper
135             # to use the perl PRNG stream by default.
136             # currently croaks because we have no binomial method
137             package Statistics::Sampler::Multinomial::DefaultPRNG {
138 2     2   8 use Carp;
  2         2  
  2         193  
139             sub new {
140 0     0     croak "No default PRNG yet implemented for Statistics::Sampler::Multinomial.\n"
141             . "Try Math::Random::MT::Auto.";
142 0           return bless {}, __PACKAGE__;
143             }
144             #sub rand {
145             # rand();
146             #}
147             sub binomial {
148             ...
149 0     0     }
150             1;
151             }
152            
153            
154             1;
155             __END__