Unbias a random generator: Difference between revisions
m
→{{header|EasyLang}}
m (→{{header|Nim}}) |
|||
(40 intermediate revisions by 23 users not shown) | |||
Line 1:
{{task}}Given a weighted one
Line 12:
This task is an implementation of [http://en.wikipedia.org/wiki/Randomness_extractor#Von_Neumann_extractor Von Neumann debiasing], first described in a 1951 paper.
<br><br>
=={{header|11l}}==
{{trans|Python}}
<syntaxhighlight lang="11l">F randN(n)
‘1,0 random generator factory with 1 appearing 1/n'th of the time’
R () -> random:(@=n) == 0
F unbiased(biased)
‘uses a biased() generator of 1 or 0, to create an unbiased one’
V (this, that) = (biased(), biased())
L this == that
(this, that) = (biased(), biased())
R this
L(n) 3..6
V biased = randN(n)
V v = (0.<1000000).map(x -> @biased())
V (v1, v0) = (v.count(1), v.count(0))
print(‘Biased(#.): count1=#., count0=#., percent=#.2’.format(n, v1, v0, 100.0 * v1 / (v1 + v0)))
v = (0.<1000000).map(x -> unbiased(@biased))
(v1, v0) = (v.count(1), v.count(0))
print(‘ Unbiased: count1=#., count0=#., percent=#.2’.format(v1, v0, 100.0 * v1 / (v1 + v0)))</syntaxhighlight>
{{out}}
<pre>
Biased(3): count1=332946, count0=667054, percent=33.29
Unbiased: count1=499751, count0=500249, percent=49.98
Biased(4): count1=250561, count0=749439, percent=25.06
Unbiased: count1=500576, count0=499424, percent=50.06
Biased(5): count1=200056, count0=799944, percent=20.01
Unbiased: count1=499975, count0=500025, percent=50.00
Biased(6): count1=165953, count0=834047, percent=16.60
Unbiased: count1=500104, count0=499896, percent=50.01
</pre>
=={{header|Ada}}==
<
procedure Bias_Unbias is
Line 64 ⟶ 100:
Ada.Text_IO.New_Line;
end loop;
end Bias_Unbias;</
Output:<pre> I Biased% UnBiased%
3 32.87 49.80
Line 74 ⟶ 110:
=={{header|Aime}}==
{{trans|C}}
<
biased(integer bias)
{
}
Line 88 ⟶ 124:
}
}
Line 99 ⟶ 135:
b = 3;
while (b <= 6) {
i = cb = cu = 0;
cb += biased(b);
cu += unbiased(b);
}
o_form("bias ~: /d2p2/%% vs /d2p2/%%\n", b, 100r * cb / n,
b += 1;
}
}</
Output:<pre>bias 3: 33.51% vs 50.27%
bias 4: 24.97% vs 49.99%
Line 124 ⟶ 156:
=={{header|AutoHotkey}}==
{{output?}}
<
Random, q, 0, 4
return q=4
Line 138 ⟶ 170:
MsgBox % "Unbiased probability of a 1 occurring: " Errorlevel/1000
StringReplace, junk, t, 1, , UseErrorLevel
MsgBox % "biased probability of a 1 occurring: " Errorlevel/1000</
=={{header|BASIC}}==
==={{header|BASIC256}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="basic256">
function randN (n)
if int(rand * n) + 1 <> 1 then return 0 else return 1
end function
function unbiased (n)
do
a = randN (n)
b = randN (n)
until a <> b
return a
end function
numveces = 100000
print "Resultados de números aleatorios sesgados e imparciales" + chr(10)
for n = 3 to 6
dim b_numveces(n) fill 0
dim u_numveces(n) fill 0
for m = 1 to numveces
x = randN (n)
b_numveces[x] += 1
x = unbiased (n)
u_numveces[x] += 1
next m
print "N = "; n
print " Biased =>", "#0="; (b_numveces[0]); " #1="; (b_numveces[1]); " ratio = "; (b_numveces[1]/numveces*100); "%"
print "Unbiased =>", "#0="; (u_numveces[0]); " #1="; (u_numveces[1]); " ratio = "; (u_numveces[1]/numveces*100); "%"
next n
end
</syntaxhighlight>
{{out}}
<pre>
Resultados de números aleatorios sesgados e imparciales
N = 3
Biased => #0=66625 #1=33375 ratio = 33.375%
Unbiased => #0=50026 #1=49974 ratio = 49.974%
N = 4
Biased => #0=74988 #1=25012 ratio = 25.012%
Unbiased => #0=49809 #1=50191 ratio = 50.191%
N = 5
Biased => #0=79893 #1=20107 ratio = 20.107%
Unbiased => #0=50102 #1=49898 ratio = 49.898%
N = 6
Biased => #0=83432 #1=16568 ratio = 16.568%
Unbiased => #0=50091 #1=49909 ratio = 49.909%
</pre>
=={{header|BBC BASIC}}==
<
biased% = 0
unbiased% = 0
Line 160 ⟶ 245:
= A%
DEF FNrandN(N%) = -(RND(N%) = 1)</
Output:
<pre>
Line 170 ⟶ 255:
=={{header|C}}==
<
#include <stdlib.h>
Line 201 ⟶ 286:
return 0;
}</
output
<pre>bias 3: 33.090% vs 49.710%
Line 207 ⟶ 292:
bias 5: 19.760% vs 49.650%
bias 6: 16.740% vs 50.030%</pre>
=={{header|C++}}==
{{trans|C#}}
<syntaxhighlight lang="cpp">#include <iostream>
#include <random>
std::default_random_engine generator;
bool biased(int n) {
std::uniform_int_distribution<int> distribution(1, n);
return distribution(generator) == 1;
}
bool unbiased(int n) {
bool flip1, flip2;
/* Flip twice, and check if the values are the same.
* If so, flip again. Otherwise, return the value of the first flip. */
do {
flip1 = biased(n);
flip2 = biased(n);
} while (flip1 == flip2);
return flip1;
}
int main() {
for (size_t n = 3; n <= 6; n++) {
int biasedZero = 0;
int biasedOne = 0;
int unbiasedZero = 0;
int unbiasedOne = 0;
for (size_t i = 0; i < 100000; i++) {
if (biased(n)) {
biasedOne++;
} else {
biasedZero++;
}
if (unbiased(n)) {
unbiasedOne++;
} else {
unbiasedZero++;
}
}
std::cout << "(N = " << n << ")\n";
std::cout << "Biased:\n";
std::cout << " 0 = " << biasedZero << "; " << biasedZero / 1000.0 << "%\n";
std::cout << " 1 = " << biasedOne << "; " << biasedOne / 1000.0 << "%\n";
std::cout << "Unbiased:\n";
std::cout << " 0 = " << unbiasedZero << "; " << unbiasedZero / 1000.0 << "%\n";
std::cout << " 1 = " << unbiasedOne << "; " << unbiasedOne / 1000.0 << "%\n";
std::cout << '\n';
}
return 0;
}</syntaxhighlight>
{{out}}
<pre>(N = 3)
Biased:
0 = 66614; 66.614%
1 = 33386; 33.386%
Unbiased:
0 = 49965; 49.965%
1 = 50035; 50.035%
(N = 4)
Biased:
0 = 75032; 75.032%
1 = 24968; 24.968%
Unbiased:
0 = 50030; 50.03%
1 = 49970; 49.97%
(N = 5)
Biased:
0 = 80178; 80.178%
1 = 19822; 19.822%
Unbiased:
0 = 49878; 49.878%
1 = 50122; 50.122%
(N = 6)
Biased:
0 = 83494; 83.494%
1 = 16506; 16.506%
Unbiased:
0 = 50085; 50.085%
1 = 49915; 49.915%</pre>
=={{header|C sharp}}==
<syntaxhighlight lang
namespace Unbias
Line 268 ⟶ 442:
}
}
}</
'''Sample Output'''
Line 287 ⟶ 461:
=={{header|Clojure}}==
<
(if (< (rand 2) (/ n)) 0 1))
Line 303 ⟶ 477:
[4 0.87684 0.5023]
[5 0.90122 0.49728]
[6 0.91526 0.5])</
=={{header|CoffeeScript}}==
<
biased_rand_function = (n) ->
# return a function that returns 0/1 with
Line 336 ⟶ 510:
stats "biased", f_biased
stats "unbiased", f_unbiased
</syntaxhighlight>
output
<pre>
Line 359 ⟶ 533:
=={{header|Common Lisp}}==
<
(defun unbiased (n)
Line 369 ⟶ 543:
(let ((u (loop repeat 10000 collect (unbiased n)))
(b (loop repeat 10000 collect (biased n))))
(format t "~a: unbiased ~d biased ~d~%" n (count 0 u) (count 0 b))))</
output
<pre>3: unbiased 4992 biased 3361
Line 377 ⟶ 551:
=={{header|D}}==
<
enum biased = (in int n) /*nothrow*/ => uniform01 < (1.0 / n);
Line 393 ⟶ 567:
M.iota.map!(_=> n.biased).sum * 100.0 / M,
M.iota.map!(_=> n.unbiased).sum * 100.0 / M);
}</
{{out}}
<pre>3: 33.441% 49.964%
Line 399 ⟶ 573:
5: 19.958% 49.987%
6: 16.660% 49.890%</pre>
=={{header|EasyLang}}==
{{trans|Java}}
<syntaxhighlight>
func biased n .
return if randomf < 1 / n
.
func unbiased n .
repeat
a = biased n
b = biased n
until a <> b
.
return a
.
m = 50000
for n = 3 to 6
c1 = 0
c2 = 0
for i to m
c1 += biased n
c2 += unbiased n
.
print n & ": " & 100 * c1 / m & " " & 100 * c2 / m
.
</syntaxhighlight>
{{out}}
<pre>
3: 33.10 50.14
4: 25.03 50.14
5: 19.81 50.04
6: 16.75 49.83
</pre>
=={{header|Elena}}==
{{trans|C#}}
ELENA 6.x :
<syntaxhighlight lang="elena">import extensions;
extension op : IntNumber
{
bool randN()
= randomGenerator.nextInt(self) == 0;
get bool Unbiased()
{
bool flip1 := self.randN();
bool flip2 := self.randN();
while (flip1 == flip2)
{
flip1 := self.randN();
flip2 := self.randN()
};
^ flip1
}
}
public program()
{
for(int n := 3; n <= 6; n += 1)
{
int biasedZero := 0;
int biasedOne := 0;
int unbiasedZero := 0;
int unbiasedOne := 0;
for(int i := 0; i < 100000; i += 1)
{
if(n.randN()) { biasedOne += 1 } else { biasedZero += 1 };
if(n.Unbiased) { unbiasedOne += 1 } else { unbiasedZero += 1 }
};
console
.printLineFormatted("(N = {0}):".padRight(17) + "# of 0"$9"# of 1"$9"% of 0"$9"% of 1", n)
.printLineFormatted("Biased:".padRight(15) + "{0}"$9"{1}"$9"{2}"$9"{3}",
biasedZero, biasedOne, biasedZero / 1000, biasedOne / 1000)
.printLineFormatted("Unbiased:".padRight(15) + "{0}"$9"{1}"$9"{2}"$9"{3}",
unbiasedZero, unbiasedOne, unbiasedZero / 1000, unbiasedOne / 1000)
}
}</syntaxhighlight>
{{out}}
<pre>
(N = 3): # of 0 # of 1 % of 0 % of 1
Biased: 66793 33207 66 33
Unbiased: 49965 50035 49 50
(N = 4): # of 0 # of 1 % of 0 % of 1
Biased: 75233 24767 75 24
Unbiased: 50106 49894 50 49
(N = 5): # of 0 # of 1 % of 0 % of 1
Biased: 80209 19791 80 19
Unbiased: 50080 49920 50 49
(N = 6): # of 0 # of 1 % of 0 % of 1
Biased: 83349 16651 83 16
Unbiased: 49699 50301 49 50
</pre>
=={{header|Elixir}}==
<
def randN(n) do
if :rand.uniform(n) == 1, do: 1, else: 0
Line 417 ⟶ 688:
ys = for _ <- 1..m, do: Random.unbiased(n)
IO.puts "#{n} #{Enum.sum(xs) / m} #{Enum.sum(ys) / m}"
end</
{{out}}
Line 429 ⟶ 700:
=={{header|ERRE}}==
<
FUNCTION RANDN(N)
Line 459 ⟶ 730:
END FOR
END PROGRAM
</syntaxhighlight>
{{out}}
<pre>N = 3 : biased = 32.66 , unbiased = 49.14
Line 468 ⟶ 739:
=={{header|Euphoria}}==
<
return rand(N) = 1
end function
Line 492 ⟶ 763:
end for
printf(1, "%d: %5.2f%% %5.2f%%\n", {b, 100 * cb / n, 100 * cu / n})
end for</
Output:
Line 502 ⟶ 773:
=={{header|F_Sharp|F#}}==
<
let random = Random()
Line 523 ⟶ 794:
printfn "%d: %5.2f%% %5.2f%%"
b (100. * float !cb / float n) (100. * float !cu / float n)
0</
{{out}}
<pre>3: 33.26% 49.97%
Line 529 ⟶ 800:
5: 19.98% 50.00%
6: 16.64% 49.69%</pre>
=={{header|Factor}}==
<syntaxhighlight lang="factor">USING: formatting kernel math math.ranges random sequences ;
IN: rosetta-code.unbias
: randN ( n -- m ) random zero? 1 0 ? ;
: unbiased ( n -- m )
dup [ randN ] dup bi 2dup = not
[ drop nip ] [ 2drop unbiased ] if ;
: test-generator ( quot -- x )
[ 1,000,000 dup ] dip replicate sum 100 * swap / ; inline
: main ( -- )
3 6 [a,b] [
dup [ randN ] [ unbiased ] bi-curry
[ test-generator ] bi@ "%d: %.2f%% %.2f%%\n" printf
] each ;
MAIN: main</syntaxhighlight>
{{out}}
<pre>
3: 33.25% 50.03%
4: 24.98% 50.02%
5: 20.03% 50.04%
6: 16.66% 49.99%
</pre>
=={{header|Fortran}}==
{{works with|Fortran|90 and later}}
<
implicit none
Line 580 ⟶ 879:
end function
end program</
Output:
<pre>3: 33.337% 49.971%
Line 586 ⟶ 885:
5: 19.971% 49.987%
6: 16.688% 50.097%</pre>
=={{header|FreeBASIC}}==
{{trans|PureBasic}}
<syntaxhighlight lang="freebasic">
Function randN (n As Ubyte) As Ubyte
If Int(Rnd * n) + 1 <> 1 Then Return 0 Else Return 1
End Function
Function unbiased (n As Ubyte) As Ubyte
Dim As Ubyte a, b
Do
a = randN (n)
b = randN (n)
Loop Until a <> b
Return a
End Function
Const count = 100000
Dim x As Ubyte
Randomize Timer
Print "Resultados de n";Chr(163);!"meros aleatorios sesgados e imparciales\n"
For n As Ubyte = 3 To 6
Dim As Integer b_count(1)
Dim As Integer u_count(1)
For m As Integer = 1 To count
x = randN (n)
b_count(x) += 1
x = unbiased (n)
u_count(x) += 1
Next m
Print "N ="; n
Print " Biased =>", "#0="; Str(b_count(0)), "#1="; Str(b_count(1)),
Print Using "ratio = ##.##%"; (b_count(1) / count * 100)
Print "Unbiased =>", "#0="; Str(u_count(0)), "#1="; Str(u_count(1)),
Print Using "ratio = ##.##%"; (u_count(1) / count * 100)
Next n
Sleep
</syntaxhighlight>
{{out}}
<pre>
Resultados de números aleatorios sesgados e imparciales
N = 3
Biased => #0=66633 #1=33367 ratio = 33.37%
Unbiased => #0=50029 #1=49971 ratio = 49.97%
N = 4
Biased => #0=75068 #1=24932 ratio = 24.93%
Unbiased => #0=50107 #1=49893 ratio = 49.89%
N = 5
Biased => #0=79998 #1=20002 ratio = 20.00%
Unbiased => #0=50049 #1=49951 ratio = 49.95%
N = 6
Biased => #0=83195 #1=16805 ratio = 16.81%
Unbiased => #0=50026 #1=49974 ratio = 49.97%
</pre>
=={{header|Fōrmulæ}}==
{{FormulaeEntry|page=https://formulae.org/?script=examples/Unbias_a_random_generator}}
'''Solution'''
[[File:Fōrmulæ - Unbias a random generator 01.png]]
[[File:Fōrmulæ - Unbias a random generator 02.png]]
[[File:Fōrmulæ - Unbias a random generator 03.png]]
'''Test cases'''
[[File:Fōrmulæ - Unbias a random generator 04.png]]
[[File:Fōrmulæ - Unbias a random generator 05.png]]
=={{header|GAP}}==
<
local v, rand;
v := [1 .. n - 1]*0;
Line 625 ⟶ 1,002:
# [ 4, 249851, 500663 ],
# [ 5, 200532, 500448 ],
# [ 6, 166746, 499859 ] ]</
=={{header|Go}}==
<
import (
Line 666 ⟶ 1,044:
u[1], u[0], float64(u[1])*100/samples)
}
}</
Output:
<pre>
Line 682 ⟶ 1,060:
=={{header|Haskell}}==
The first task:
<
import Control.Monad
import Text.Printf
randN :: MonadRandom m => Int -> m Int
randN n = fromList [(0, fromIntegral n-1), (1, 1)]</
Examples of use:
Line 696 ⟶ 1,074:
The second task. Returns the unbiased generator for any given random generator.
<
unbiased g = do x <- g
y <- g
if x /= y then return y else unbiased g</
Examples of use:
Line 708 ⟶ 1,086:
The third task:
<
where
showCounts b = do
Line 715 ⟶ 1,093:
printf "n = %d biased: %d%% unbiased: %d%%\n" b r1 r2
counts g = (`div` 100) . length . filter (== 1) <$> replicateM 10000 g</
Output:
Line 729 ⟶ 1,107:
This solution works in both languages. Both <tt>randN</tt> and
<tt>unbiased</tt> are generators in the Icon/Unicon sense.
<
iters := \A[1] | 10000
write("ratios of 0 to 1 from ",iters," trials:")
Line 758 ⟶ 1,136:
procedure randN(n)
repeat suspend if 1 = ?n then 1 else 0
end</
and a sample run:
<pre>->ubrn 100000
Line 773 ⟶ 1,151:
=={{header|J}}==
<
unbiased=: i.@# { ::$: 2 | 0 3 -.~ _2 #.\ 4&* randN@# ]</
Example use:
<
1 0 0 0 1 0 0 0 0 0
unbiased 10#6
1 0 0 1 0 0 1 0 1 1</
Some example counts (these are counts of the number of 1s which appear in a test involving 100 random numbers):
<
30
+/randN 100#4
Line 800 ⟶ 1,178:
49
+/unbiased 100#6
47</
Note that these results are
=={{header|Java}}==
<
public static boolean biased(int n) {
return Math.random() < 1.0 / n;
Line 831 ⟶ 1,209:
}
}
}</
Output:
<pre>3: 33,11% 50,23%
Line 837 ⟶ 1,215:
5: 20.05% 50.00%
6: 17.00% 49.88%</pre>
=={{header|jq}}==
'''Adapted from [[#Wren]]'''
{{works with|jq}}
'''Also works with gojq, the Go implementation of jq'''
In this entry, /dev/random is used as a source of entropy,
and so a suitable invocation would be:
<pre>
< /dev/random tr -cd '0-9' | fold -w 1 | jq -Mcnr -f unbias.jq
</pre>
'''unbias.jq'''
<syntaxhighlight lang=jq>
### Utility Functions
# Output: a PRN in range(0;$n) where $n is .
def prn:
if . == 1 then 0
else . as $n
| (($n-1)|tostring|length) as $w
| [limit($w; inputs)] | join("") | tonumber
| if . < $n then . else ($n | prn) end
end;
def round: ((. * 100) | floor) / 100;
# input: n
# output: boolean, such that P(true) == 1/n
def biased:
prn == 1 | debug;
def unbiased:
. as $n
| {}
| until( .a != .b; {a: ($n|biased), b: ($n|biased)})
| .a;
def task($m):
def f(x;y;z): "\(x): \(y|round)% \(z|round)%";
range(3;7) as $n
| reduce range(0; $m) as $i ( {c1:0, c2:0};
if ($n|biased) then .c1 += 1 else . end
| if ($n|unbiased) then .c2 += 1 else . end)
| f($n; 100 * .c1 / $m; 100 * .c2 / $m);
task(50000)
</syntaxhighlight>
{{output}}
<pre>
3: 33.61% 50.09%
4: 25.48% 50.21%
5: 20.32% 50.47%
6: 16.35% 49.72%
</pre>
=={{header|Julia}}==
<syntaxhighlight lang="julia">using Printf
randN(N) = () -> rand(1:N) == 1 ? 1 : 0
function unbiased(biased::Function)
this, that = biased(), biased()
while this == that this, that = biased(), biased() end
return this
end
@printf "%2s | %10s | %5s | %5s | %8s" "N" "bias./unb." "1s" "0s" "pct ratio"
const nrep = 10000
for N in 3:6
biased = randN(N)
v = collect(biased() for __ in 1:nrep)
v1, v0 = count(v .== 1), count(v .== 0)
@printf("%2i | %10s | %5i | %5i | %5.2f%%\n", N, "biased", v1, v0, 100 * v1 / nrep)
v = collect(unbiased(biased) for __ in 1:nrep)
v1, v0 = count(v .== 1), count(v .== 0)
@printf("%2i | %10s | %5i | %5i | %5.2f%%\n", N, "unbiased", v1, v0, 100 * v1 / nrep)
end</syntaxhighlight>
{{out}}
<pre> N | bias./unb. | 1s | 0s | pct ratio
3 | biased | 3286 | 6714 | 32.86%
3 | unbiased | 4986 | 5014 | 49.86%
4 | biased | 2473 | 7527 | 24.73%
4 | unbiased | 4986 | 5014 | 49.86%
5 | biased | 1992 | 8008 | 19.92%
5 | unbiased | 5121 | 4879 | 51.21%
6 | biased | 1663 | 8337 | 16.63%
6 | unbiased | 5040 | 4960 | 50.40%</pre>
=={{header|Kotlin}}==
{{trans|Java}}
<
fun biased(n: Int) = Math.random() < 1.0 / n
Line 867 ⟶ 1,335:
println(f.format(n, 100.0 * c1 / m, 100.0 * c2 / m))
}
}</
Sample output:
Line 875 ⟶ 1,343:
5: 19.91% 50.07%
6: 16.71% 50.14%
</pre>
=={{header|Lua}}==
<
local function randN(n)
return function()
Line 990 ⟶ 1,389:
demonstrate(100000)
</syntaxhighlight>
Output:
Line 1,008 ⟶ 1,407:
</pre>
=={{header|Mathematica}}/{{header|Wolfram Language}}==
<
unbiased[bias_, n_] := DeleteCases[rand[bias, {n, 2}], {a_, a_}][[All, 1]]
count = 1000000;
TableForm[
Table[{n, Total[rand[n, count]]/count // N,
Total[#]/Length[#] &@unbiased[n, count] // N}, {n, 3, 6}],
TableHeadings -> {None, {n, "biased", "unbiased"}}]</syntaxhighlight>
{{out}}
<pre>n biased unbiased
3 0.33312 0.500074
4 0.24932 0.499883
5 0.1998 0.498421
6 0.16620 0.49805</pre>
=={{header|NetRexx}}==
{{trans|Java}}
<
options replace format comments java crossref symbols binary
Line 1,065 ⟶ 1,459:
end n
return
</syntaxhighlight>
'''Output:'''
<pre>
Line 1,076 ⟶ 1,470:
=={{header|Nim}}==
{{trans|Python}}
<
type randProc = proc: range[0..1]
randomize()
proc randN(n
result = proc: range[0..1] = ord(
proc unbiased(biased: randProc): range[0..1] =
result = biased()
var that = biased()
Line 1,099 ⟶ 1,491:
var cnt0, cnt1 = 0
for x in v:
if x == 0: inc cnt0 else: inc cnt1
echo &"Biased({n}) → count1 = {cnt1}, count0 = {cnt0}, percent = {100*cnt1 / (cnt1+cnt0):.3f}"
v = newSeqWith(1_000_000, unbiased(biased))
Line 1,108 ⟶ 1,498:
cnt1 = 0
for x in v:
if x == 0: inc cnt0 else: inc cnt1
echo &"Unbiased → count1 = {cnt1}, count0 = {cnt0}, percent = {100*cnt1 / (cnt1+cnt0):.3f}"</syntaxhighlight>
{{out}}
<pre>Biased(3) → count1 = 249654, count0 = 750346, percent = 24.965
Biased(4) →
Biased(5) →
Biased(6) →
=={{header|OCaml}}==
<
if Random.int n = 0 then 1 else 0
Line 1,143 ⟶ 1,532:
Printf.printf "%d: %5.2f%% %5.2f%%\n"
b (100.0 *. float !cb /. float n) (100.0 *. float !cu /. float n)
done</
Output:
Line 1,154 ⟶ 1,543:
=={{header|PARI/GP}}==
GP's random number generation is high-quality, using Brent's [http://maths.anu.edu.au/~brent/random.html XORGEN]. Thus this program is slow: the required 400,000 unbiased numbers generated through this bias/unbias scheme take nearly a second. This requires about two million calls to <code>random</code>, which in turn generate a total of about three million calls to the underlying random number generator through the rejection strategy. The overall efficiency of the scheme is 0.8% for 32-bit and 0.4% for 64-bit...
<
unbiased(N)={
my(a,b);
Line 1,163 ⟶ 1,552:
)
};
for(n=3,6,print(n"\t"sum(k=1,1e5,unbiased(n))"\t"sum(k=1,1e5,randN(n))))</
Output:
Line 1,172 ⟶ 1,561:
=={{header|Perl}}==
<
my $n = shift;
return int(rand($n) / ($n - 1));
Line 1,192 ⟶ 1,581:
100 * sqrt($fixed[0] * $fixed[1]) / ($fixed[0] + $fixed[1])**1.5);
}</
Output:
<pre>Bias 3: 6684 3316, 66.84+-0.471% fixed: 2188 2228, 49.5471+-0.752%
Line 1,198 ⟶ 1,587:
Bias 5: 7993 2007, 79.93+-0.401% fixed: 1564 1597, 49.478+-0.889%
Bias 6: 8309 1691, 83.09+-0.375% fixed: 1403 1410, 49.8756+-0.943%</pre>
=={{header|Phix}}==
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">randN</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">N</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">return</span> <span style="color: #7060A8;">rand</span><span style="color: #0000FF;">(</span><span style="color: #000000;">N</span><span style="color: #0000FF;">)</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">unbiased</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">N</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">while</span> <span style="color: #004600;">true</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">a</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">randN</span><span style="color: #0000FF;">(</span><span style="color: #000000;">N</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">!=</span><span style="color: #000000;">randN</span><span style="color: #0000FF;">(</span><span style="color: #000000;">N</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">a</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">n</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">100000</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">=</span><span style="color: #000000;">3</span> <span style="color: #008080;">to</span> <span style="color: #000000;">6</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">cb</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">cu</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">n</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">cb</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">randN</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">cu</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">unbiased</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"%d: biased:%.0f%% unbiased:%.0f%%\n"</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">b</span><span style="color: #0000FF;">,(</span><span style="color: #000000;">cb</span><span style="color: #0000FF;">/</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">100</span><span style="color: #0000FF;">,(</span><span style="color: #000000;">cu</span><span style="color: #0000FF;">/</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">100</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<!--</syntaxhighlight>-->
{{out}}
<small>''(Rounded to the nearest whole percent, of course)''</small>
<pre>
3:
4:
5: biased:20
6:
</pre>
=={{header|PicoLisp}}==
<
(if (= 1 (rand 1 N)) 1 0) )
Line 1,272 ⟶ 1,634:
(setq A (randN N))
(setq B (randN N)) ) )
A ) )</
Test:
<
(tab (2 1 7 2 7 2)
N ":"
Line 1,284 ⟶ 1,646:
(let S 0 (do 10000 (inc 'S (unbiased N))))
2 )
"%" ) )</
Output:
<pre> 3: 33.21 % 50.48 %
Line 1,292 ⟶ 1,654:
=={{header|PL/I}}==
<syntaxhighlight lang="pl/i">
test: procedure options (main); /* 20 Nov. 2012 */
Line 1,322 ⟶ 1,684:
end test;
</syntaxhighlight>
Results:
<pre>
Line 1,334 ⟶ 1,696:
=={{header|PowerShell}}==
{{works with|PowerShell|2}}
<syntaxhighlight lang="powershell">
function randN ( [int]$N )
{
Line 1,350 ⟶ 1,712:
return $X
}
</syntaxhighlight>
Note: The [pscustomobject] type accelerator, used to simplify making the test output look pretty, requires version 3.0 or higher.
<syntaxhighlight lang="powershell">
$Tests = 1000
ForEach ( $N in 3..6 )
Line 1,368 ⟶ 1,730:
"Unbiased Ones out of $Test" = $Unbiased }
}
</syntaxhighlight>
{{out}}
<pre>
Line 1,380 ⟶ 1,742:
=={{header|PureBasic}}==
<
If Random(n) <> 1
ProcedureReturn 0
Line 1,414 ⟶ 1,776:
output + #tab$ + " ratio=" + StrF(u_count(1) / #count * 100, 2) + "%" + #LF$
Next
MessageRequester("Biased and Unbiased random number results", output)</
Sample output:
<pre>---------------------------
Line 1,433 ⟶ 1,795:
=={{header|Python}}==
<
import random
Line 1,460 ⟶ 1,822:
v = [unbiased(biased) for x in range(1000000)]
v1, v0 = v.count(1), v.count(0)
print ( " Unbiased = %r" % (Stats(v1, v0, 100. * v1/(v1 + v0)), ) )</
'''Sample output'''
Line 1,472 ⟶ 1,834:
Biased(6) = Stats(count1=167561, count0=832439, percent=16.7561)
Unbiased = Stats(count1=499963, count0=500037, percent=49.996299999999998)</pre>
=={{header|Quackery}}==
<syntaxhighlight lang="quackery"> $ "bigrat.qky" loadfile
[ random 0 = ] is randN ( n --> n )
[ dup randN
over randN
2dup = iff
2drop again
drop nip ] is unbias ( n --> n )
[ dup echo say " biased --> "
0
1000000 times
[ over randN if 1+ ]
nip 1000000 6 point$ echo$ ] is showbias ( n --> )
[ dup echo say " unbiased --> "
0
1000000 times
[ over unbias if 1+ ]
nip 1000000 6 point$ echo$ ] is showunbias ( n --> )
' [ 3 4 5 6 ]
witheach
[ dup cr
showbias cr
showunbias cr ] </syntaxhighlight>
{{out}}
<pre>3 biased --> 0.333225
3 unbiased --> 0.500147
4 biased --> 0.249658
4 unbiased --> 0.499851
5 biased --> 0.200169
5 unbiased --> 0.500073
6 biased --> 0.166804
6 unbiased --> 0.499045
</pre>
=={{header|R}}==
<
unbiased = function(f)
Line 1,485 ⟶ 1,892:
N = N,
biased = mean(replicate(samples, randN(N))),
unbiased = mean(replicate(samples, unbiased(function() randN(N)))))))))</
Sample output:
Line 1,496 ⟶ 1,903:
=={{header|Racket}}==
<
#lang racket
;; Using boolean #t/#f instead of 1/0
Line 1,510 ⟶ 1,917:
(printf "Count: ~a => Biased: ~a%; Unbiased: ~a%.\n"
n (try% biased) (try% (unbiased biased))))
</syntaxhighlight>
{{out}}
<pre>
Line 1,518 ⟶ 1,925:
Count: 6 => Biased: 17%; Unbiased: 50%.
</pre>
=={{header|Raku}}==
(formerly Perl 6)
{{trans|Perl}}
{{works with|Rakudo|2020.08.1}}
<syntaxhighlight lang="raku" line>sub randN ( $n where 3..6 ) {
return ( $n.rand / ($n - 1) ).Int;
}
sub unbiased ( $n where 3..6 ) {
my $n1;
repeat { $n1 = randN($n) } until $n1 != randN($n);
return $n1;
}
my $iterations = 1000;
for 3 .. 6 -> $n {
my ( @raw, @fixed );
for ^$iterations {
@raw[ randN($n) ]++;
@fixed[ unbiased($n) ]++;
}
printf "N=%d randN: %s, %4.1f%% unbiased: %s, %4.1f%%\n",
$n, map { .raku, .[1] * 100 / $iterations }, @raw, @fixed;
}</syntaxhighlight>
Output:<pre>N=3 randN: [676, 324], 32.4% unbiased: [517, 483], 48.3%
N=4 randN: [734, 266], 26.6% unbiased: [489, 511], 51.1%
N=5 randN: [792, 208], 20.8% unbiased: [494, 506], 50.6%
N=6 randN: [834, 166], 16.6% unbiased: [514, 486], 48.6%</pre>
=={{header|REXX}}==
<
parse arg # R seed . /*
if #=='' | #=="," then #=1000 /*#
if R=='' | R=="," then R=6 /*R
if datatype(seed, 'W') then call random ,,seed /*
say left('',5) ctr("N",5) ctr(@b) ctr(@b'%') ctr(@ub) ctr(@ub"%") ctr('samples')
dash=
do N=3 to R;
do j=1 for #; b=b +
say
end /*N*/
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
pct:
randN:
unbiased:
<pre>
──N── ───biased─── ──biased%─── ──unbiased── ─unbiased%── ──samples───
Line 1,548 ⟶ 1,984:
6 178 17.80% 488 48.80% 1000
</pre>
<pre>
──N── ───biased─── ──biased%─── ──unbiased── ─unbiased%── ──samples───
Line 1,556 ⟶ 1,992:
6 1644 16.44% 4982 49.82% 10000
</pre>
<pre>
──N── ───biased─── ──biased%─── ──unbiased── ─unbiased%── ──samples───
Line 1,590 ⟶ 2,026:
=={{header|Ring}}==
<
for n = 3 to 6
biased = 0
Line 1,610 ⟶ 2,046:
m = (random(m) = 1)
return m
</syntaxhighlight>
Output:
<pre>
Line 1,617 ⟶ 2,053:
N = 5 : biased = 16.65%, unbiased = 48.86%
N = 6 : biased = 13.31%, unbiased = 49.96%
</pre>
=={{header|RPL}}==
≪ INV RAND ≥ ≫ ‘'''RandN'''’ STO
≪ 0 DUP '''WHILE''' DUP2 == '''REPEAT'''
DROP2 DUP '''RandN''' OVER '''RandN'''
'''END''' ROT DROP2
≫ ‘'''Unbiased'''’ STO
≪ 3 6 '''FOR''' n
(0,0) 1 10000 '''START'''
n '''RandN''' + n '''Unbiased''' i * →NUM +
'''NEXT''' 10000 / '''NEXT'''
≫ ‘SHOW’ STO
{{in}}
<pre>
SHOW
</pre>
{{out}}
<pre>
4: (0.3272,0.4976)
3: (0.2504,0.4961)
2: (0.201,0.5008)
1: (0.166,0.4952)
</pre>
=={{header|Ruby}}==
<
rand(bias) == 0 ? 1 : 0
end
Line 1,641 ⟶ 2,102:
counter[:bias] = bias
puts counter.values_at(*keys).join("\t")
end</
{{output}}
<pre>
Line 1,649 ⟶ 2,110:
5 199767 500354
6 166163 499809
</pre>
=={{header|Rust}}==
<syntaxhighlight lang="rust">#![feature(inclusive_range_syntax)]
extern crate rand;
use rand::Rng;
fn rand_n<R: Rng>(rng: &mut R, n: u32) -> usize {
rng.gen_weighted_bool(n) as usize // maps `false` to 0 and `true` to 1
}
fn unbiased<R: Rng>(rng: &mut R, n: u32) -> usize {
let mut bit = rand_n(rng, n);
while bit == rand_n(rng, n) {
bit = rand_n(rng, n);
}
bit
}
fn main() {
const SAMPLES: usize = 100_000;
let mut rng = rand::weak_rng();
println!(" Bias rand_n unbiased");
for n in 3..=6 {
let mut count_biased = 0;
let mut count_unbiased = 0;
for _ in 0..SAMPLES {
count_biased += rand_n(&mut rng, n);
count_unbiased += unbiased(&mut rng, n);
}
let b_percentage = 100.0 * count_biased as f64 / SAMPLES as f64;
let ub_percentage = 100.0 * count_unbiased as f64 / SAMPLES as f64;
println!(
"bias {}: {:0.2}% {:0.2}%",
n, b_percentage, ub_percentage
);
}
}</syntaxhighlight>
{{output}}
<pre>
Bias rand_n unbiased
bias 3: 33.32% 49.80%
bias 4: 25.22% 50.16%
bias 5: 19.91% 50.00%
bias 6: 16.66% 49.95%
</pre>
=={{header|Scala}}==
<
def unbiased( n:Int ) = { def loop : Boolean = { val a = biased(n); if( a != biased(n) ) a else loop }; loop }
Line 1,664 ⟶ 2,174:
"%d: %2.2f%% %2.2f%%".format(i, 100.0*c1/m, 100.0*c2/m)
}</
{{output}}
<pre>3: 33.09% 49.79%
Line 1,672 ⟶ 2,182:
=={{header|Seed7}}==
<
include "float.s7i";
Line 1,705 ⟶ 2,215:
" " <& flt(100 * sumUnbiased) / flt(tests) digits 3 lpad 6);
end for;
end func;</
Output:
Line 1,716 ⟶ 2,226:
=={{header|Sidef}}==
{{trans|
<
n.rand / (n-1) -> int
}
Line 1,738 ⟶ 2,248:
printf("N=%d randN: %s, %4.1f%% unbiased: %s, %4.1f%%\n",
n, [raw, fixed].map {|a| (a.dump, a[1] * 100 / iterations) }...)
}</
{{out}}
<pre>
Line 1,748 ⟶ 2,258:
=={{header|Tcl}}==
<
proc randN n {expr {rand()*$n < 1}}
Line 1,770 ⟶ 2,280:
puts [format "unbiased %d => #0=%d #1=%d ratio=%.2f%%" $n $c(0) $c(1) \
[expr {100.*$c(1)/$i}]]
}</
Sample output:
<pre>
Line 1,781 ⟶ 2,291:
biased 6 => #0=833623 #1=166377 ratio=16.64%
unbiased 6 => #0=500518 #1=499482 ratio=49.95%
</pre>
=={{header|Visual Basic .NET}}==
{{trans|C#}}
<syntaxhighlight lang="vbnet">Module Module1
Dim random As New Random()
Function RandN(n As Integer) As Boolean
Return random.Next(0, n) = 0
End Function
Function Unbiased(n As Integer) As Boolean
Dim flip1 As Boolean
Dim flip2 As Boolean
Do
flip1 = RandN(n)
flip2 = RandN(n)
Loop While flip1 = flip2
Return flip1
End Function
Sub Main()
For n = 3 To 6
Dim biasedZero = 0
Dim biasedOne = 0
Dim unbiasedZero = 0
Dim unbiasedOne = 0
For i = 1 To 100000
If RandN(n) Then
biasedOne += 1
Else
biasedZero += 1
End If
If Unbiased(n) Then
unbiasedOne += 1
Else
unbiasedZero += 1
End If
Next
Console.WriteLine("(N = {0}):".PadRight(17) + "# of 0" + vbTab + "# of 1" + vbTab + "% of 0" + vbTab + "% of 1", n)
Console.WriteLine("(Biased: {0}):".PadRight(15) + "{0}" + vbTab + "{1}" + vbTab + "{2}" + vbTab + "{3}", biasedZero, biasedOne, biasedZero / 1000, biasedOne / 1000)
Console.WriteLine("(UnBiased: {0}):".PadRight(15) + "{0}" + vbTab + "{1}" + vbTab + "{2}" + vbTab + "{3}", unbiasedZero, unbiasedOne, unbiasedZero / 1000, unbiasedOne / 1000)
Next
End Sub
End Module</syntaxhighlight>
{{out}}
<pre>(N = 3): # of 0 # of 1 % of 0 % of 1
(Biased: 66844): 66844 33156 66.844 33.156
(UnBiased: 50081):50081 49919 50.081 49.919
(N = 4): # of 0 # of 1 % of 0 % of 1
(Biased: 74827): 74827 25173 74.827 25.173
(UnBiased: 50095):50095 49905 50.095 49.905
(N = 5): # of 0 # of 1 % of 0 % of 1
(Biased: 79878): 79878 20122 79.878 20.122
(UnBiased: 50544):50544 49456 50.544 49.456
(N = 6): # of 0 # of 1 % of 0 % of 1
(Biased: 83449): 83449 16551 83.449 16.551
(UnBiased: 49907):49907 50093 49.907 50.093</pre>
=={{header|Wren}}==
{{trans|Kotlin}}
{{libheader|Wren-fmt}}
<syntaxhighlight lang="wren">import "random" for Random
import "./fmt" for Fmt
var rand = Random.new()
var biased = Fn.new { |n| rand.float() < 1 / n }
var unbiased = Fn.new { |n|
while (true) {
var a = biased.call(n)
var b = biased.call(n)
if (a != b) return a
}
}
var m = 50000
var f = "$d: $2.2f\% $2.2f\%"
for (n in 3..6) {
var c1 = 0
var c2 = 0
for (i in 0...m) {
if (biased.call(n)) c1 = c1 + 1
if (unbiased.call(n)) c2 = c2 + 1
}
Fmt.print(f, n, 100 * c1 / m, 100 * c2 / m)
}</syntaxhighlight>
{{out}}
Sample output:
<pre>
3: 33.36% 49.91%
4: 24.80% 49.89%
5: 20.07% 49.95%
6: 16.64% 50.17%
</pre>
=={{header|zkl}}==
<
fcn unbiased(randN){ while((a:=randN())==randN()){} a }</
<
foreach N in ([3..6]){
"%d: biased: %3.2f%%, unbiased: %3.2f%%".fmt(N,
Line 1,792 ⟶ 2,402:
(0).reduce(Z,'wrap(s,_){ s+unbiased(randN.fp(N)) },0.0)/Z*100)
.println();
}</
{{out}}
<pre>
|