Unbias a random generator: Difference between revisions
m
→{{header|EasyLang}}
(added Lua example) |
|||
| (138 intermediate revisions by 67 users not shown) | |||
Line 1:
{{task}}Given a weighted one
;Task details:
* Use your language's random number generator to create a function/method/subroutine/... '''randN''' that returns a one or a zero, but with one occurring, on average, 1 out of N times, where N is an integer from the range 3 to 6 inclusive.
* Create a function '''unbiased''' that uses only randN as its source of randomness to become an unbiased generator of random ones and zeroes.
* For N over its range, generate and show counts of the outputs of randN and unbiased(randN).
<br>
The actual unbiasing should be done by generating two numbers at a time from randN and only returning a 1 or 0 if they are different. As long as you always return the first number or always return the second number, the probabilities discussed above should take over the biased probability of randN.
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}}==
<syntaxhighlight lang="ada">with Ada.Text_IO; with Ada.Numerics.Discrete_Random;
procedure Bias_Unbias is
Modulus: constant Integer := 60; -- lcm of {3,4,5,6}
type M is mod Modulus;
package Rand is new Ada.Numerics.Discrete_Random(M);
Gen: Rand.Generator;
subtype Bit is Integer range 0 .. 1;
function Biased_Bit(Bias_Base: Integer) return Bit is
begin
if (Integer(Rand.Random(Gen))* Bias_Base) / Modulus > 0 then
return 0;
else
return 1;
end if;
end Biased_Bit;
function Unbiased_Bit(Bias_Base: Integer) return Bit is
A, B: Bit := 0;
begin
while A = B loop
A := Biased_Bit(Bias_Base);
B := Biased_Bit(Bias_Base);
end loop;
return A;
end Unbiased_Bit;
package FIO is new Ada.Text_IO.Float_IO(Float);
Counter_B, Counter_U: Natural;
Number_Of_Samples: constant Natural := 10_000;
begin
Rand.Reset(Gen);
Ada.Text_IO.Put_Line(" I Biased% UnBiased%");
for I in 3 .. 6 loop
Counter_B := 0;
Counter_U := 0;
for J in 1 .. Number_Of_Samples loop
Counter_B := Counter_B + Biased_Bit(I);
Counter_U := Counter_U + Unbiased_Bit(I);
end loop;
Ada.Text_IO.Put(Integer'Image(I));
FIO.Put(100.0 * Float(Counter_B) / Float(Number_Of_Samples), 5, 2, 0);
FIO.Put(100.0 * Float(Counter_U) / Float(Number_Of_Samples), 5, 2, 0);
Ada.Text_IO.New_Line;
end loop;
end Bias_Unbias;</syntaxhighlight>
Output:<pre> I Biased% UnBiased%
3 32.87 49.80
4 24.49 50.22
5 19.73 50.05
6 16.75 50.19
</pre>
=={{header|Aime}}==
{{trans|C}}
<syntaxhighlight lang="aime">integer
biased(integer bias)
{
1 ^ min(drand(bias - 1), 1);
}
integer
unbiased(integer bias)
{
integer a;
while ((a = biased(bias)) == biased(bias)) {
}
a;
}
integer
main(void)
{
integer b, n, cb, cu, i;
n = 10000;
b = 3;
while (b <= 6) {
i = cb = cu = 0;
while ((i += 1) <= n) {
cb += biased(b);
cu += unbiased(b);
}
o_form("bias ~: /d2p2/%% vs /d2p2/%%\n", b, 100r * cb / n,
100r * cu / n);
b += 1;
}
0;
}</syntaxhighlight>
Output:<pre>bias 3: 33.51% vs 50.27%
bias 4: 24.97% vs 49.99%
bias 5: 19.93% vs 49.92%
bias 6: 16.32% vs 49.36%</pre>
=={{header|AutoHotkey}}==
{{output?}}
<syntaxhighlight lang="autohotkey">Biased(){
Random, q, 0, 4
return q=4
}
Unbiased(){
Loop
If ((a := Biased()) != biased())
return a
}
Loop 1000
t .= biased(), t2 .= unbiased()
StringReplace, junk, t2, 1, , UseErrorLevel
MsgBox % "Unbiased probability of a 1 occurring: " Errorlevel/1000
StringReplace, junk, t, 1, , UseErrorLevel
MsgBox % "biased probability of a 1 occurring: " Errorlevel/1000</syntaxhighlight>
=={{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}}==
<syntaxhighlight lang="bbcbasic"> FOR N% = 3 TO 6
biased% = 0
unbiased% = 0
FOR I% = 1 TO 10000
IF FNrandN(N%) biased% += 1
IF FNunbiased(N%) unbiased% += 1
NEXT
PRINT "N = ";N% " : biased = "; biased%/100 "%, unbiased = "; unbiased%/100 "%"
NEXT
END
DEF FNunbiased(N%)
LOCAL A%,B%
REPEAT
A% = FNrandN(N%)
B% = FNrandN(N%)
UNTIL A%<>B%
= A%
DEF FNrandN(N%) = -(RND(N%) = 1)</syntaxhighlight>
Output:
<pre>
N = 3 : biased = 33.57%, unbiased = 49.94%
N = 4 : biased = 25.34%, unbiased = 50.76%
N = 5 : biased = 20.06%, unbiased = 50.04%
N = 6 : biased = 16.25%, unbiased = 50.13%
</pre>
=={{header|C}}==
<syntaxhighlight lang="c">#include <stdio.h>
#include <stdlib.h>
int biased(int bias)
{
/* balance out the bins, being pedantic */
int r, rand_max = RAND_MAX - (RAND_MAX % bias);
while ((r = rand()) > rand_max);
return r < rand_max / bias;
}
int unbiased(int bias)
{
int a;
while ((a = biased(bias)) == biased(bias));
return a;
}
int main()
{
int b, n = 10000, cb, cu, i;
for (b = 3; b <= 6; b++) {
for (i = cb = cu = 0; i < n; i++) {
cb += biased(b);
cu += unbiased(b);
}
printf("bias %d: %5.3f%% vs %5.3f%%\n", b,
100. * cb / n, 100. * cu / n);
}
return 0;
}</syntaxhighlight>
output
<pre>bias 3: 33.090% vs 49.710%
bias 4: 25.130% vs 49.430%
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="csharp">using System;
namespace Unbias
{
internal class Program
{
private static void Main(string[] args)
{
// Demonstrate.
for (int n = 3; n <= 6; n++)
{
int biasedZero = 0, biasedOne = 0, unbiasedZero = 0, unbiasedOne = 0;
for (int i = 0; i < 100000; i++)
{
if (randN(n))
biasedOne++;
else
biasedZero++;
if (Unbiased(n))
unbiasedOne++;
else
unbiasedZero++;
}
Console.WriteLine("(N = {0}):".PadRight(17) + "# of 0\t# of 1\t% of 0\t% of 1", n);
Console.WriteLine("Biased:".PadRight(15) + "{0}\t{1}\t{2}\t{3}",
biasedZero, biasedOne,
biasedZero/1000, biasedOne/1000);
Console.WriteLine("Unbiased:".PadRight(15) + "{0}\t{1}\t{2}\t{3}",
unbiasedZero, unbiasedOne,
unbiasedZero/1000, unbiasedOne/1000);
}
}
private static 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 = randN(n);
flip2 = randN(n);
} while (flip1 == flip2);
return flip1;
}
private static readonly Random random = new Random();
private static bool randN(int n)
{
// Has an 1/n chance of returning 1. Otherwise it returns 0.
return random.Next(0, n) == 0;
}
}
}</syntaxhighlight>
'''Sample Output'''
<pre>(N = 3): # of 0 # of 1 % of 0 % of 1
Biased: 66867 33133 66 33
Unbiased: 49843 50157 49 50
(N = 4): # of 0 # of 1 % of 0 % of 1
Biased: 74942 25058 74 25
Unbiased: 50192 49808 50 49
(N = 5): # of 0 # of 1 % of 0 % of 1
Biased: 80203 19797 80 19
Unbiased: 49928 50072 49 50
(N = 6): # of 0 # of 1 % of 0 % of 1
Biased: 83205 16795 83 16
Unbiased: 49744 50256 49 50
</pre>
=={{header|Clojure}}==
<
(if (< (rand 2) (/ n)) 0 1))
Line 24 ⟶ 477:
[4 0.87684 0.5023]
[5 0.90122 0.49728]
[6 0.91526 0.5])</
=={{header|CoffeeScript}}==
<syntaxhighlight lang="coffeescript">
biased_rand_function = (n) ->
# return a function that returns 0/1 with
# 1 appearing only 1/Nth of the time
cap = 1/n
->
if Math.random() < cap
1
else
0
unbiased_function = (f) ->
->
while true
[n1, n2] = [f(), f()]
return n1 if n1 + n2 == 1
stats = (label, f) ->
cnt = 0
sample_size = 10000000
for i in [1...sample_size]
cnt += 1 if f() == 1
console.log "ratio of 1s: #{cnt / sample_size} [#{label}]"
for n in [3..6]
console.log "\n---------- n = #{n}"
f_biased = biased_rand_function(n)
f_unbiased = unbiased_function f_biased
stats "biased", f_biased
stats "unbiased", f_unbiased
</syntaxhighlight>
output
<pre>
> coffee unbiased.coffee
---------- n = 3
ratio of 1s: 0.3333343 [biased]
ratio of 1s: 0.4999514 [unbiased]
---------- n = 4
ratio of 1s: 0.2499751 [biased]
ratio of 1s: 0.4998067 [unbiased]
---------- n = 5
ratio of 1s: 0.199729 [biased]
ratio of 1s: 0.5003183 [unbiased]
---------- n = 6
ratio of 1s: 0.1664843 [biased]
ratio of 1s: 0.4997813 [unbiased]
</pre>
=={{header|Common Lisp}}==
<syntaxhighlight lang="lisp">(defun biased (n) (if (zerop (random n)) 0 1))
(defun unbiased (n)
(loop with x do
(if (/= (setf x (biased n)) (biased n))
(return x))))
(loop for n from 3 to 6 do
(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))))</syntaxhighlight>
output
<pre>3: unbiased 4992 biased 3361
4: unbiased 4988 biased 2472
5: unbiased 5019 biased 1987
6: unbiased 4913 biased 1658</pre>
=={{header|D}}==
<
while ((a = bias.biased) == bias.biased) {}
return a;
}
void main() {
foreach (immutable n; 3 .. 7)
}</syntaxhighlight>
{{out}}
<pre>3: 33.441% 49.964%
4: 24.953% 49.910%
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}}==
<syntaxhighlight lang="elixir">defmodule Random do
def randN(n) do
if :rand.uniform(n) == 1, do: 1, else: 0
end
def unbiased(n) do
{x, y} = {randN(n), randN(n)}
if x != y, do: x, else: unbiased(n)
end
end
IO.puts "N biased unbiased"
m = 10000
for n <- 3..6 do
xs = for _ <- 1..m, do: Random.randN(n)
ys = for _ <- 1..m, do: Random.unbiased(n)
IO.puts "#{n} #{Enum.sum(xs) / m} #{Enum.sum(ys) / m}"
end</syntaxhighlight>
{{out}}
<pre>
N biased unbiased
3 0.3356 0.5043
4 0.2523 0.4996
5 0.2027 0.5041
6 0.1647 0.4912
</pre>
=={{header|ERRE}}==
<syntaxhighlight lang="erre">PROGRAM UNBIAS
FUNCTION RANDN(N)
RANDN=INT(1+N*RND(1))=1
END FUNCTION
PROCEDURE UNBIASED(N->RIS)
LOCAL A,B
REPEAT
A=RANDN(N)
B=RANDN(N)
UNTIL A<>B
RIS=A
END PROCEDURE
BEGIN
PRINT(CHR$(12);) ! CLS
RANDOMIZE(TIMER)
FOR N=3 TO 6 DO
BIASED=0
UNBIASED=0
FOR I=1 TO 10000 DO
IF RANDN(N) THEN biased+=1
UNBIASED(N->RIS)
IF RIS THEN unbiased+=+1
END FOR
PRINT("N =";N;" : biased =";biased/100;", unbiased =";unbiased/100)
END FOR
END PROGRAM
</syntaxhighlight>
{{out}}
<pre>N = 3 : biased = 32.66 , unbiased = 49.14
N = 4 : biased = 25.49 , unbiased = 49.92
N = 5 : biased = 20.53 , unbiased = 50
N = 6 : biased = 17.43 , unbiased = 50.43
</pre>
=={{header|Euphoria}}==
<syntaxhighlight lang="euphoria">function randN(integer N)
return rand(N) = 1
end function
function unbiased(integer N)
integer a
while 1 do
a = randN(N)
if a != randN(N) then
return a
end if
end while
end function
constant n = 10000
integer cb, cu
for b = 3 to 6 do
cb = 0
cu = 0
for i = 1 to n do
cb += randN(b)
cu += unbiased(b)
end for
printf(1, "%d: %5.2f%% %5.2f%%\n", {b, 100 * cb / n, 100 * cu / n})
end for</syntaxhighlight>
Output:
<pre>3: 33.
4:
5:
6: 16.
</pre>
=={{header|F_Sharp|F#}}==
<syntaxhighlight lang="fsharp">open System
let random = Random()
let randN = random.Next >> (=)0 >> Convert.ToInt32
let rec unbiased n =
let a = randN n
if a <> randN n then a else unbiased n
[<EntryPoint>]
let main argv =
let n = if argv.Length > 0 then UInt32.Parse(argv.[0]) |> int else 100000
for b = 3 to 6 do
let cb = ref 0
let cu = ref 0
for i = 1 to n do
cb := !cb + randN b
cu := !cu + unbiased b
printfn "%d: %5.2f%% %5.2f%%"
b (100. * float !cb / float n) (100. * float !cu / float n)
0</syntaxhighlight>
{{out}}
<pre>3: 33.26% 49.97%
4: 25.02% 50.22%
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 109 ⟶ 879:
end function
end program</
Output:
<pre>3: 33.337% 49.971%
Line 115 ⟶ 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}}==
<syntaxhighlight lang="gap">RandNGen := function(n)
local v, rand;
v := [1 .. n - 1]*0;
Add(v, 1);
rand := function()
return Random(v);
end;
return rand;
end;
UnbiasedGen := function(rand)
local unbiased;
unbiased := function()
local a, b;
while true do
a := rand();
b := rand();
if a <> b then
break;
fi;
od;
return a;
end;
return unbiased;
end;
range := [2 .. 6];
v := List(range, RandNGen);
w := List(v, UnbiasedGen);
apply := gen -> Sum([1 .. 1000000], n -> gen());
# Some tests (2 is added as a witness, since in this case RandN is already unbiased)
PrintArray(TransposedMat([range, List(v, apply), List(w, apply)]));
# [ [ 2, 499991, 499041 ],
# [ 3, 333310, 500044 ],
# [ 4, 249851, 500663 ],
# [ 5, 200532, 500448 ],
# [ 6, 166746, 499859 ] ]</syntaxhighlight>
=={{header|Go}}==
<syntaxhighlight lang="go">package main
import (
"fmt"
"math/rand"
)
const samples = 1e6
func main() {
fmt.Println("Generator 1 count 0 count % 1 count")
for n := 3; n <= 6; n++ {
// function randN, per task description
randN := func() int {
if rand.Intn(n) == 0 {
return 1
}
return 0
}
var b [2]int
for x := 0; x < samples; x++ {
b[randN()]++
}
fmt.Printf("randN(%d) %7d %7d %5.2f%%\n",
n, b[1], b[0], float64(b[1])*100/samples)
// function unbiased, per task description
unbiased := func() (b int) {
for b = randN(); b == randN(); b = randN() {
}
return
}
var u [2]int
for x := 0; x < samples; x++ {
u[unbiased()]++
}
fmt.Printf("unbiased %7d %7d %5.2f%%\n",
u[1], u[0], float64(u[1])*100/samples)
}
}</syntaxhighlight>
Output:
<pre>
Generator 1 count 0 count % 1 count
randN(3) 332711 667289 33.27%
unbiased 499649 500351 49.96%
randN(4) 249742 750258 24.97%
unbiased 499434 500566 49.94%
randN(5) 200318 799682 20.03%
unbiased 499100 500900 49.91%
randN(6) 166900 833100 16.69%
unbiased 499973 500027 50.00%
</pre>
=={{header|Haskell}}==
The first task:
<syntaxhighlight lang="haskell">import Control.Monad.Random
import Control.Monad
import Text.Printf
randN :: MonadRandom m => Int -> m Int
randN n = fromList [(0, fromIntegral n-1), (1, 1)]</syntaxhighlight>
Examples of use:
<pre>λ> replicateM 20 (randN 2)
[0,0,1,0,0,1,0,1,1,0,0,1,1,1,1,0,1,1,0,0]
λ> replicateM 20 (randN 5)
[0,1,0,0,0,0,0,0,1,0,0,0,0,1,0,0,1,0,1,0]</pre>
The second task. Returns the unbiased generator for any given random generator.
<syntaxhighlight lang="haskell">unbiased :: (MonadRandom m, Eq x) => m x -> m x
unbiased g = do x <- g
y <- g
if x /= y then return y else unbiased g</syntaxhighlight>
Examples of use:
<pre>λ> replicateM 20 (unbiased (randN 5))
[0,0,1,0,1,1,1,0,0,0,1,1,1,0,1,1,0,0,1,0]
λ> replicateM 20 (unbiased (fromList [(True,10),(False,1)]))
[True,True,False,True,True,True,False,True,False,True,True,False,False,True,False,True,True,False,False,True]</pre>
The third task:
<syntaxhighlight lang="haskell">main = forM_ [3..6] showCounts
where
showCounts b = do
r1 <- counts (randN b)
r2 <- counts (unbiased (randN b))
printf "n = %d biased: %d%% unbiased: %d%%\n" b r1 r2
counts g = (`div` 100) . length . filter (== 1) <$> replicateM 10000 g</syntaxhighlight>
Output:
<pre>
n = 3 biased: 33% unbiased: 49%
n = 4 biased: 24% unbiased: 50%
n = 5 biased: 19% unbiased: 50%
n = 6 biased: 16% unbiased: 49%
</pre>
=={{header|Icon}} and {{header|Unicon}}==
This solution works in both languages. Both <tt>randN</tt> and
<tt>unbiased</tt> are generators in the Icon/Unicon sense.
<syntaxhighlight lang="unicon">procedure main(A)
iters := \A[1] | 10000
write("ratios of 0 to 1 from ",iters," trials:")
every n := 3 to 6 do {
results_randN := table(0)
results_unbiased := table(0)
every 1 to iters do {
results_randN[randN(n)] +:= 1
results_unbiased[unbiased(n)] +:= 1
}
showResults(n, "randN", results_randN)
showResults(n, "unbiased", results_unbiased)
}
end
procedure showResults(n, s, t)
write(n," ",left(s,9),":",t[0],"/",t[1]," = ",t[0]/real(t[1]))
end
procedure unbiased(n)
repeat {
n1 := randN(n)
n2 := randN(n)
if n1 ~= n2 then suspend n1
}
end
procedure randN(n)
repeat suspend if 1 = ?n then 1 else 0
end</syntaxhighlight>
and a sample run:
<pre>->ubrn 100000
ratios of 0 to 1 from 100000 trials:
3 randN :66804/33196 = 2.012411133871551
3 unbiased :49812/50188 = 0.9925081692834941
4 randN :75017/24983 = 3.002721850858584
4 unbiased :50000/50000 = 1.0
5 randN :79990/20010 = 3.997501249375312
5 unbiased :50073/49927 = 1.002924269433373
6 randN :83305/16695 = 4.989817310572027
6 unbiased :49911/50089 = 0.9964463255405378
-></pre>
=={{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 144 ⟶ 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 175 ⟶ 1,209:
}
}
}</
Output:
<pre>3: 33,11% 50,23%
Line 181 ⟶ 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}}
<syntaxhighlight lang="scala">// version 1.1.2
fun biased(n: Int) = Math.random() < 1.0 / n
fun unbiased(n: Int): Boolean {
var a: Boolean
var b: Boolean
do {
a = biased(n)
b = biased(n)
}
while (a == b)
return a
}
fun main(args: Array<String>) {
val m = 50_000
val f = "%d: %2.2f%% %2.2f%%"
for (n in 3..6) {
var c1 = 0
var c2 = 0
for (i in 0 until m) {
if (biased(n)) c1++
if (unbiased(n)) c2++
}
println(f.format(n, 100.0 * c1 / m, 100.0 * c2 / m))
}
}</syntaxhighlight>
Sample output:
<pre>
3: 33.19% 50.19%
4: 25.29% 49.85%
5: 19.91% 50.07%
6: 16.71% 50.14%
</pre>
=={{header|Lua}}==
<
local function randN(n)
return function()
if math.random() <
end
end
Line 228 ⟶ 1,389:
demonstrate(100000)
</syntaxhighlight>
Output:
Line 244 ⟶ 1,405:
biased 83407 16593 83.407 16.593
unbias 49820 50180 49.82 50.18
</pre>
=={{header|Mathematica}}/{{header|Wolfram Language}}==
<syntaxhighlight lang="mathematica">rand[bias_, n_] := 1 - Unitize@RandomInteger[bias - 1, n]
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}}
<syntaxhighlight lang="netrexx">/* NetRexx */
options replace format comments java crossref symbols binary
runSample(arg)
return
-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
method biased(n = int) public static returns boolean
return Math.random() < 1.0 / n
-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
method unbiased(n = int) public static returns boolean
a = boolean
b = boolean
loop until a \= b
a = biased(n)
b = biased(n)
end
return a
-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
method runSample(arg) private static
parse arg Mx .
if Mx.length <= 0 then Mx = 50000
M = int Mx
loop n = int 3 to 6
c1 = int 0
c2 = int 0
loop for M
if biased(n) then c1 = c1 + 1
if unbiased(n) then c2 = c2 + 1
end
say Rexx(n).right(3)':' Rexx(100.0 * c1 / M).format(6, 2)'%' Rexx(100.0 * c2 / M).format(6, 2)'%'
end n
return
</syntaxhighlight>
'''Output:'''
<pre>
3: 32.78% 49.98%
4: 24.72% 50.31%
5: 19.95% 50.34%
6: 17.20% 50.20%
</pre>
=={{header|Nim}}==
{{trans|Python}}
<syntaxhighlight lang="nim">import random, sequtils, strformat
type randProc = proc: range[0..1]
randomize()
proc randN(n: Positive): randProc =
result = proc: range[0..1] = ord(rand(n) == 0)
proc unbiased(biased: randProc): range[0..1] =
result = biased()
var that = biased()
while result == that:
result = biased()
that = biased()
for n in 3..6:
var biased = randN(n)
var v = newSeqWith(1_000_000, biased())
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))
cnt0 = 0
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
Unbiased → count1 = 498631, count0 = 501369, percent = 49.863
Biased(4) → count1 = 200117, count0 = 799883, percent = 20.012
Unbiased → count1 = 499425, count0 = 500575, percent = 49.943
Biased(5) → count1 = 167218, count0 = 832782, percent = 16.722
Unbiased → count1 = 500020, count0 = 499980, percent = 50.002
Biased(6) → count1 = 143388, count0 = 856612, percent = 14.339
Unbiased → count1 = 500378, count0 = 499622, percent = 50.038</pre>
=={{header|OCaml}}==
<syntaxhighlight lang="ocaml">let randN n =
if Random.int n = 0 then 1 else 0
let rec unbiased n =
let a = randN n in
if a <> randN n then a else unbiased n
let () =
Random.self_init();
let n = 50_000 in
for b = 3 to 6 do
let cb = ref 0 in
let cu = ref 0 in
for i = 1 to n do
cb := !cb + (randN b);
cu := !cu + (unbiased b);
done;
Printf.printf "%d: %5.2f%% %5.2f%%\n"
b (100.0 *. float !cb /. float n) (100.0 *. float !cu /. float n)
done</syntaxhighlight>
Output:
<pre>3: 33.07% 49.90%
4: 25.11% 49.85%
5: 19.82% 50.09%
6: 16.51% 50.51%</pre>
=={{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...
<syntaxhighlight lang="parigp">randN(N)=!random(N);
unbiased(N)={
my(a,b);
while(1,
a=randN(N);
b=randN(N);
if(a!=b, return(a))
)
};
for(n=3,6,print(n"\t"sum(k=1,1e5,unbiased(n))"\t"sum(k=1,1e5,randN(n))))</syntaxhighlight>
Output:
<pre>3 49997 33540
4 49988 24714
5 50143 20057
6 49913 16770</pre>
=={{header|Perl}}==
<syntaxhighlight lang="perl">sub randn {
my $n = shift;
return int(rand($n) / ($n - 1));
}
for my $n (3 .. 6) {
print "Bias $n: ";
my (@raw, @fixed);
for (1 .. 10000) {
my $x = randn($n);
$raw[$x]++;
$fixed[$x]++ if randn($n) != $x
}
print "@raw, ";
printf("%3g+-%.3g%%\tfixed: ", $raw[0]/100,
100 * sqrt($raw[0] * $raw[1]) / ($raw[0] + $raw[1])**1.5);
print "@fixed, ";
printf("%3g+-%.3g%%\n", 100*$fixed[0]/($fixed[0] + $fixed[1]),
100 * sqrt($fixed[0] * $fixed[1]) / ($fixed[0] + $fixed[1])**1.5);
}</syntaxhighlight>
Output:
<pre>Bias 3: 6684 3316, 66.84+-0.471% fixed: 2188 2228, 49.5471+-0.752%
Bias 4: 7537 2463, 75.37+-0.431% fixed: 1924 1845, 51.048+-0.814%
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: biased:33% unbiased:50%
4: biased:25% unbiased:50%
5: biased:20% unbiased:50%
6: biased:17% unbiased:50%
</pre>
=={{header|PicoLisp}}==
<
(if (= 1 (rand 1 N)) 1 0) )
Line 256 ⟶ 1,634:
(setq A (randN N))
(setq B (randN N)) ) )
A ) )</
Test:
<
(tab (2 1 7 2 7 2)
N ":"
Line 268 ⟶ 1,646:
(let S 0 (do 10000 (inc 'S (unbiased N))))
2 )
"%" ) )</
Output:
<pre> 3: 33.21 % 50.48 %
Line 274 ⟶ 1,652:
5: 20.04 % 49.75 %
6: 16.32 % 49.02 %</pre>
=={{header|PL/I}}==
<syntaxhighlight lang="pl/i">
test: procedure options (main); /* 20 Nov. 2012 */
randN: procedure(N) returns (bit (1));
declare N fixed (1);
declare random builtin;
declare r fixed (2) external initial (-1);
if r >= 0 then do; r = r-1; return ('0'b); end;
r = random()*2*N;
return ('1'b);
end randN;
random: procedure returns (bit(1));
declare (r1, r2) bit (1);
do until (r1 ^= r2);
r1 = randN(N); r2 = randN(N);
end;
return (r1);
end random;
declare (biasedrn, unbiasedrn) (100) bit (1);
declare N fixed (1);
put ('N Biased Unbiased (tally of 100 random numbers)');
do N = 3 to 6;
biasedrn = randN(N); unbiasedrn = random;
put skip edit (N, sum(biasedrn), sum(unbiasedrn)) (F(1), 2 F(10));
end;
end test;
</syntaxhighlight>
Results:
<pre>
N Biased Unbiased (tally of 100 random numbers)
3 24 42
4 18 47
5 16 41
6 11 53
</pre>
=={{header|PowerShell}}==
{{works with|PowerShell|2}}
<syntaxhighlight lang="powershell">
function randN ( [int]$N )
{
[int]( ( Get-Random -Maximum $N ) -eq 0 )
}
function unbiased ( [int]$N )
{
do {
$X = randN $N
$Y = randN $N
}
While ( $X -eq $Y )
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 )
{
$Biased = 0
$Unbiased = 0
ForEach ( $Test in 1..$Tests )
{
$Biased += randN $N
$Unbiased += unbiased $N
}
[pscustomobject]@{ N = $N
"Biased Ones out of $Test" = $Biased
"Unbiased Ones out of $Test" = $Unbiased }
}
</syntaxhighlight>
{{out}}
<pre>
N Biased Ones out of 1000 Unbiased Ones out of 1000
- ----------------------- -------------------------
3 322 503
4 273 518
5 217 515
6 173 486
</pre>
=={{header|PureBasic}}==
<
If Random(n) <> 1
ProcedureReturn 0
Line 309 ⟶ 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 328 ⟶ 1,795:
=={{header|Python}}==
<
import random
def randN(N):
" 1,0 random generator factory with 1 appearing 1/N'th of the time"
return lambda: random.randrange(N) == 0
def unbiased(biased):
Line 359 ⟶ 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 371 ⟶ 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}}==
<syntaxhighlight lang="rsplus">randN = function(N) sample.int(N, 1) == 1
unbiased = function(f)
{while ((x <- f()) == f()) {}
x}
samples = 10000
print(t(round(d = 2, sapply(3:6, function(N) c(
N = N,
biased = mean(replicate(samples, randN(N))),
unbiased = mean(replicate(samples, unbiased(function() randN(N)))))))))</syntaxhighlight>
Sample output:
<pre> N biased unbiased
[1,] 3 0.32 0.50
[2,] 4 0.24 0.50
[3,] 5 0.21 0.49
[4,] 6 0.16 0.51</pre>
=={{header|Racket}}==
<syntaxhighlight lang="racket">
#lang racket
;; Using boolean #t/#f instead of 1/0
(define ((randN n)) (zero? (random n)))
(define ((unbiased biased))
(let loop () (let ([r (biased)]) (if (eq? r (biased)) (loop) r))))
;; Counts
(define N 1000000)
(for ([n (in-range 3 7)])
(define (try% R) (round (/ (for/sum ([i N]) (if (R) 1 0)) N 1/100)))
(define biased (randN n))
(printf "Count: ~a => Biased: ~a%; Unbiased: ~a%.\n"
n (try% biased) (try% (unbiased biased))))
</syntaxhighlight>
{{out}}
<pre>
Count: 3 => Biased: 33%; Unbiased: 50%.
Count: 4 => Biased: 25%; Unbiased: 50%.
Count: 5 => Biased: 20%; Unbiased: 50%.
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}}==
<syntaxhighlight lang="rexx">/*REXX program generates unbiased random numbers and displays the results to terminal.*/
parse arg # R seed . /*obtain optional arguments from the CL*/
if #=='' | #=="," then #=1000 /*#: the number of SAMPLES to be used.*/
if R=='' | R=="," then R=6 /*R: the high number for the range. */
if datatype(seed, 'W') then call random ,,seed /*Specified? Then use for RANDOM seed.*/
dash='─'; @b="biased"; @ub='un'@b /*literals for the SAY column headers. */
say left('',5) ctr("N",5) ctr(@b) ctr(@b'%') ctr(@ub) ctr(@ub"%") ctr('samples')
dash=
do N=3 to R; b=0; u=0
do j=1 for #; b=b + randN(N); u=u + unbiased()
end /*j*/
say left('', 5) ctr(N, 5) ctr(b) pct(b) ctr(u) pct(u) ctr(#)
end /*N*/
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
ctr: return center( arg(1), word(arg(2) 12, 1), left(dash, 1)) /*show hdr│numbers.*/
pct: return ctr( format(arg(1) / # * 100, , 2)'%' ) /*2 decimal digits.*/
randN: parse arg z; return random(1, z)==z /*ret 1 if rand==Z.*/
unbiased: do until x\==randN(N); x=randN(N); end; return x /* " unbiased RAND*/</syntaxhighlight>
{{out|output|text= when using the default inputs:}}
<pre>
──N── ───biased─── ──biased%─── ──unbiased── ─unbiased%── ──samples───
3 348 34.80% 541 54.10% 1000
4 259 25.90% 479 47.90% 1000
5 188 18.80% 475 47.50% 1000
6 178 17.80% 488 48.80% 1000
</pre>
{{out|output|text= when using the input of: <tt> 10000 </tt>}}
<pre>
──N── ───biased─── ──biased%─── ──unbiased── ─unbiased%── ──samples───
3 3435 34.35% 4995 49.95% 10000
4 2535 25.35% 4957 49.57% 10000
5 2019 20.19% 4958 49.58% 10000
6 1644 16.44% 4982 49.82% 10000
</pre>
{{out|output|text= when using the input of: <tt> 100000 30 </tt>}}
<pre>
──N── ───biased─── ──biased%─── ──unbiased── ─unbiased%── ──samples───
3 33301 33.30% 50066 50.07% 100000
4 25359 25.36% 49401 49.40% 100000
5 20026 20.03% 49966 49.97% 100000
6 16579 16.58% 49956 49.96% 100000
7 14294 14.29% 50008 50.01% 100000
8 12402 12.40% 50479 50.48% 100000
9 11138 11.14% 50099 50.10% 100000
10 9973 9.97% 49988 49.99% 100000
11 9062 9.06% 50009 50.01% 100000
12 8270 8.27% 49929 49.93% 100000
13 7704 7.70% 49876 49.88% 100000
14 7223 7.22% 50414 50.41% 100000
15 6725 6.73% 50043 50.04% 100000
16 6348 6.35% 50252 50.25% 100000
17 5900 5.90% 49977 49.98% 100000
18 5583 5.58% 49991 49.99% 100000
19 5139 5.14% 49958 49.96% 100000
20 4913 4.91% 50198 50.20% 100000
21 4714 4.71% 49892 49.89% 100000
22 4517 4.52% 49760 49.76% 100000
23 4226 4.23% 50021 50.02% 100000
24 4174 4.17% 50141 50.14% 100000
25 4005 4.01% 49816 49.82% 100000
26 3890 3.89% 49819 49.82% 100000
27 3705 3.71% 50036 50.04% 100000
28 3567 3.57% 49665 49.67% 100000
29 3481 3.48% 50094 50.09% 100000
30 3355 3.36% 49831 49.83% 100000
</pre>
=={{header|Ring}}==
<syntaxhighlight lang="ring">
for n = 3 to 6
biased = 0
unb = 0
for i = 1 to 10000
biased += randN(n)
unb += unbiased(n)
next
see "N = " + n + " : biased = " + biased/100 + "%, unbiased = " + unb/100 + "%" + nl
next
func unbiased nr
while 1
a = randN(nr)
if a != randN(nr) return a ok
end
func randN m
m = (random(m) = 1)
return m
</syntaxhighlight>
Output:
<pre>
N = 3 : biased = 25.38%, unbiased = 50.12%
N = 4 : biased = 20.34%, unbiased = 49.17%
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}}==
<syntaxhighlight lang="ruby">def rand_n(bias)
rand(bias) == 0 ? 1 : 0
end
def unbiased(bias)
a, b = rand_n(bias), rand_n(bias) until a != b #loop until a and b are 0,1 or 1,0
a
end
runs = 1_000_000
keys = %i(bias biased unbiased) #use [:bias,:biased,:unbiased] in Ruby < 2.0
puts keys.join("\t")
(3..6).each do |bias|
counter = Hash.new(0) # counter will respond with 0 when key is not known
runs.times do
counter[:biased] += 1 if rand_n(bias) == 1 #the first time, counter has no key for :biased, so it will respond 0
counter[:unbiased] += 1 if unbiased(bias) == 1
end
counter[:bias] = bias
puts counter.values_at(*keys).join("\t")
end</syntaxhighlight>
{{output}}
<pre>
bias biased unbiased
3 333043 500161
4 249133 499393
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}}==
<syntaxhighlight lang="scala">def biased( n:Int ) = scala.util.Random.nextFloat < 1.0 / n
def unbiased( n:Int ) = { def loop : Boolean = { val a = biased(n); if( a != biased(n) ) a else loop }; loop }
for( i <- (3 until 7) ) println {
val m = 50000
var c1,c2 = 0
(0 until m) foreach { j => if( biased(i) ) c1 += 1; if( unbiased(i) ) c2 += 1 }
"%d: %2.2f%% %2.2f%%".format(i, 100.0*c1/m, 100.0*c2/m)
}</syntaxhighlight>
{{output}}
<pre>3: 33.09% 49.79%
4: 24.92% 49.92%
5: 19.75% 49.92%
6: 16.67% 50.23%</pre>
=={{header|Seed7}}==
<syntaxhighlight lang="seed7">$ include "seed7_05.s7i";
include "float.s7i";
const func integer: randN (in integer: n) is
return ord(rand(1, n) = 1);
const func integer: unbiased (in integer: n) is func
result
var integer: unbiased is 0;
begin
repeat
unbiased := randN(n);
until unbiased <> randN(n);
end func;
const proc: main is func
local
const integer: tests is 50000;
var integer: n is 0;
var integer: sumBiased is 0;
var integer: sumUnbiased is 0;
var integer: count is 0;
begin
for n range 3 to 6 do
sumBiased := 0;
sumUnbiased := 0;
for count range 1 to tests do
sumBiased +:= randN(n);
sumUnbiased +:= unbiased(n);
end for;
writeln(n <& ": " <& flt(100 * sumBiased) / flt(tests) digits 3 lpad 6 <&
" " <& flt(100 * sumUnbiased) / flt(tests) digits 3 lpad 6);
end for;
end func;</syntaxhighlight>
Output:
<pre>
3: 33.004 50.024
4: 25.158 50.278
5: 20.186 49.978
6: 16.570 49.936
</pre>
=={{header|Sidef}}==
{{trans|Raku}}
<syntaxhighlight lang="ruby">func randN (n) {
n.rand / (n-1) -> int
}
func unbiased(n) {
var n1 = nil
do { n1 = randN(n) } while (n1 == randN(n))
return n1
}
var iterations = 1000
for n in (3..6) {
var raw = []
var fixed = []
iterations.times {
raw[ randN(n) ] := 0 ++
fixed[ unbiased(n) ] := 0 ++
}
printf("N=%d randN: %s, %4.1f%% unbiased: %s, %4.1f%%\n",
n, [raw, fixed].map {|a| (a.dump, a[1] * 100 / iterations) }...)
}</syntaxhighlight>
{{out}}
<pre>
N=3 randN: [661, 339], 33.9% unbiased: [497, 503], 50.3%
N=4 randN: [765, 235], 23.5% unbiased: [493, 507], 50.7%
N=5 randN: [812, 188], 18.8% unbiased: [509, 491], 49.1%
N=6 randN: [820, 180], 18.0% unbiased: [510, 490], 49.0%
</pre>
=={{header|Tcl}}==
<
proc randN n {expr {rand()*$n < 1}}
Line 395 ⟶ 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 407 ⟶ 2,292:
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}}==
<syntaxhighlight lang="zkl">fcn randN(N){ (not (0).random(N)).toInt() }
fcn unbiased(randN){ while((a:=randN())==randN()){} a }</syntaxhighlight>
<syntaxhighlight lang="zkl">const Z=0d100_000;
foreach N in ([3..6]){
"%d: biased: %3.2f%%, unbiased: %3.2f%%".fmt(N,
(0).reduce(Z,'wrap(s,_){ s+randN(N) },0.0)/Z*100,
(0).reduce(Z,'wrap(s,_){ s+unbiased(randN.fp(N)) },0.0)/Z*100)
.println();
}</syntaxhighlight>
{{out}}
<pre>
3: biased: 33.46%, unbiased: 49.80%
4: biased: 24.95%, unbiased: 50.01%
5: biased: 19.89%, unbiased: 50.18%
6: biased: 16.75%, unbiased: 50.22%
</pre>
{{omit from|GUISS}}
[[Category:Randomness]]
| |||