みたいになるかと思います。このおみくじは5種類だからそれほどではないけど、これが100種類とか1000000種類とかだとどうでしょう。種類数に比例した数の比較演算が必要で時間がかかります。O(N)ですね。(ちょっと工夫して二分探索を使うと O(logN) になります。)$r = rand(1) if ($r < 0.15) { return "大吉"; } elsif ($r < 0.45) { return "中吉"; } elsif ($r < 0.75) { return "吉"; } elsif ($r < 0.95) { return "凶"; } else { return "大凶"; }
#!/usr/bin/perl use strict; use warnings; ###### 別名法 (alias method) ### 取りうる値の個数が有限個の任意の離散分布 ### P(X=x_{k}) = p_{k} (1 <= k <= K) ### に従う乱数を発生させるためにウォーカー(Walker, A.J.)が考案した ### 巧妙な方法である。逆関数法と違って、K がいくら大きくても ### 1個の乱数の発生に必要な比較演算が1回だけであるという特徴を持つ。 ### ただし、そのために最初に次のように表を用意しておく必要がある。 ### (参考文献:"自然科学の統計学", 東京大学出版会, 1992) my @ns = @ARGV; my $K = @ns; # K my $sum = 0; foreach my $i (@ns) { $sum += $i; } my @p = map {$_ / $sum} @ns; # p_{k} print "p_{k} :@p\n\n"; ### [初期設定] 次の手順に従って a_{k}, v_{k} (1 <= k <= K) を求める ### 1 : v_{k} ← K*p_{k} (1 <= k <= K) my @v = map {$K * $_} @p; # v_{k} ### 2 : v_{k} >= 1 を満たす k の集合を G, ### v_{k} < 1 を満たす k の集合を S とする。 my @G; my @S; foreach my $k (0..$K-1) { push @G, $k if $v[$k] >= 1; push @S, $k if $v[$k] < 1; } ### 3 : S が空集合でない限り、3.1 から 3.4 までを繰り返し実行する。 my @a; # a_{k} while (@S and @G) { ### 3.1 : S, G から1つずつ任意の要素を選ぶ。それを k, l とする。 my $k = $S[0]; my $l = $G[0]; print "v_{k}: @v\n"; print "S: @S\n"; print "G: @G\n"; print "k=$k, l=$l\n"; print "\n"; ### 3.2 : a_{k} ← x_{l}, v_{l} ← v_{l} - (1 - v_{k}) $a[$k] = $l; $v[$l] = $v[$l] - (1 - $v[$k]); ### 3.3 : S から k を除く。 shift @S; ### 3.4 : v_{l} < 1 なら l を G から S に移す。 if ($v[$l] < 1) { push @S, (shift @G); } } print "k p_{k} Kp_{k} v_{k} a_{k}\n"; for (my $i = 0; $i < $K; $i++) { print "$i $p[$i] ".($K*$p[$i])." $v[$i] ".(defined $a[$i]?$a[$i]:"*")."\n"; } print "\n"; ### [乱数発生] sub gen { my ($a_r, $v_r) = @_; ### 1 : (0, K) 上の一様乱数 V = KU を発生する。 my $V = rand(@$v_r); ### 2 : k ← L V 」+ 1, u ← k - V my $k = int($V); my $u = $V - $k; ### 3 : u <= v_{k} ならば X ← x_{k}, そうでなければ X ← a_{k} if ($u <= $v_r->[$k]) { return $k; } else { return $a_r->[$k]; } } ### 実験 my $N = 1000000; my @count; for (my $i = 0; $i < $N; $i++) { my $X = gen(\@a, \@v); $count[$X]++; } for (my $i = 0; $i < $K; $i++) { print "$i ".("*"x($count[$i]/10000))." $count[$i]\n"; }
% aliasmethod.pl 0.12 0.24 0.38 0.16 0.10 p_{k} :0.12 0.24 0.38 0.16 0.1 v_{k}: 0.6 1.2 1.9 0.8 0.5 S: 0 3 4 G: 1 2 k=0, l=1 v_{k}: 0.6 0.8 1.9 0.8 0.5 S: 3 4 1 G: 2 k=3, l=2 v_{k}: 0.6 0.8 1.7 0.8 0.5 S: 4 1 G: 2 k=4, l=2 v_{k}: 0.6 0.8 1.2 0.8 0.5 S: 1 G: 2 k=1, l=2 k p_{k} Kp_{k} v_{k} a_{k} 0 0.12 0.6 0.6 1 1 0.24 1.2 0.8 2 2 0.38 1.9 1 * 3 0.16 0.8 0.8 2 4 0.1 0.5 0.5 2 0 *********** 119632 1 ************************ 240086 2 ************************************* 379950 3 *************** 159991 4 ********** 100341
#!/usr/bin/perl use strict; use warnings; my @ns = @ARGV; my ($a_r, $v_r) = pre(@ns+0, \@ns); my $N = 1000000; my @c = do {my @t; map {$t[gen($a_r, $v_r)]++} 1..$N; @t}; print join("", map {"$_ ".("*"x($c[$_]/($N/100)))." $c[$_]\n"} 0..$#ns); sub pre { my ($K, $ns_r) = @_; my $sum = do {my $s; map {$s += $_} @$ns_r; $s}; my (@v, @a, @p, @G, @S); map {$_ /= $sum; push @p, $_; push @v, $K * $_} @$ns_r; map {$v[$_] >= 1 ? push @G, $_ : push @S, $_} 0..$K-1; for (my ($k, $l); $k = shift @S, $l = $a[$k] = $G[0]; ) { $v[$l] = $v[$l] - (1 - $v[$k]); push @S, shift @G if $v[$l] < 1; } return \@a, \@v; } sub gen { my ($a_r, $v_r) = @_; my $V = rand(@$v_r); my $k = int($V); return $V - $k <= $v_r->[$k] ? $k : $a_r->[$k]; }
別名法(alias method)。「取りうる値の個数(n), が有限個の離散分布に従う乱数を, nに無関係なO(1)の速度で生成できる, 簡単で効率的な方法である. ただし, 初期設定にO(n)の時間とO(n)のメモリを必要とする.」
In computing, the alias method is a family of efficient algorithms for sampling from a discrete probability distribution, due to A. J. Walker. The algorithms typically use O(n log n) or O(n) preprocessing time, after which random values can be drawn from the distribution in O(1) time.