Random numbers

[combike]

Background: I'm trying to study something in probabilities, and to do so, I need a lot of random numbers in the range of 1 - 3.
The problem: The RND() function doesn't seem all that random.
What I tried: I used RANDOMIZE to initialize, then SAMPLE = INT(RND(1) * 3) + 1  to generate the random numbers.
What I found: What I got were consistently more 1s than 2s and more 2s than 3s.  Not by a huge amount, but there was a definite trend.  The breakdown was always around 33.6%. 33.3%, and 33.0% (roughly - I know that adds up to 99.9%).  It doesn't seem significant, but it was very consistent.  In a group of 1 million numbers, it's quite noticeable.
My solution: I decided to use the first digit after the decimal point of the random number MOD 3 to get the values I needed, and because the MOD operation will result in an unbalance of results, I eliminate all the 9s.
The code:
    do
      initial = rnd(1)
      temp = initial * 10
      sample = int(temp)
      sample = (temp mod 3) + 1
   loop while initial >= 0.9
This results in a more even distribution.
The question: Does anyone see anything wrong with what I've done, or is the RND() function kind of skewed?

[B+]

I've noticed the same in the very start of using RND(seed) when 0 or 1, it seems the Random Number Generator starts off on low side. It gets better as you go. I've even tried a similar fix.

Have you tried running off many random numbers and then taking distribution check? how much is many?



seed = .7
for i = 1 to 100000
   r = rnd(seed)
next


for i= 1 to 100000
   r = int(rnd(seed)*10) + 1
   d(r) = d(r) + 1
next
for i = 1 to 10
   print d(i), 100*d(i)/100000;"%"
   IF i < 6 THEN low = low + d(i) ELSE high = high + d(i)
next
PRINT
PRINT 100 * low / 100000, 100 * high / 100000




Append:
There I tried, the lower numbers still seem favored over the higher.

Modified code with seed since first post.


edit 2018-05-08: Attention! Of RND(X) the (X) part is meaningless! It is not a seed for a consistent set of "random" numbers. Use Randomize for that.
Ha! No wonder the "running off" a few numbers first did not work!

[Rod]

If you want random numbers don’t use randomise. Randomise initiates a perfect sequence of similar numbers. Some folks seem to think randomise “randomises” a sequence, it does not it just starts the random sequence at the same place. So if it appears skewed it will always appear skewed.

[tsh73]

Yes JB RND() function is known to be biased.

In cases that it does matter:
1) there is a trick to get unbiased 0/1 from biased generator - you can use these bits to make number you need
2) there are known formulas of RNG renerators, for example
    Microsoft Linear congruential generator
(this one is easy and fast)

Links:
justbasic.wikispaces.com/BB___RND#RND()-Example-RND bias and "unbiased coin" function
rosettacode.org/wiki/Unbias_a_random_generator
rosettacode.org/wiki/Linear_congruential_generator#Liberty_BASIC

Some code
N=30000'0
dim a(3)

print "Build-in RND (skewed)"
for i= 1 TO N
SCAN
   sample= INT(RND(1) * 3) + 1
   a(sample)=a(sample)+1
next
print a(1), a(2), a(3)

print "Combike's code"
redim a(3)  'clear
for i= 1 TO N
SCAN
  do
     initial = rnd(1)
     temp = initial * 10
     sample = int(temp)
     'sample = (temp mod 3) + 1    'wrong
     sample = (sample  mod 3) + 1
  loop while initial >= 0.9
   a(sample )=a(sample )+1
next
print a(1), a(2), a(3)

print "Unbiased coin - 2 coins, use as bits"
redim a(3)  'clear
for i= 1 TO N
SCAN
   do
       do
           n1 = rnd(1)>0.5
           n2 = rnd(1)>0.5
       loop while n1=n2
       coin1 = n1
       do
           n1 = rnd(1)>0.5
           n2 = rnd(1)>0.5
       loop while n1=n2
       coin2 = n1
       sample = coin1+coin2*2+1    'coin1, coun2 is unbiased coins. Use them as bits
   loop while sample>3             're-roll if got out of 1,2,3

   a(sample )=a(sample )+1
next
print a(1), a(2), a(3)

print "Unbiased coin - function" 'same, looks cleaner, but calling functions costs time
redim a(3)  'clear
for i= 1 TO N
SCAN
   do
       coin1 = unBiasedCoin()
       coin2 = unBiasedCoin()
       sample = coin1+coin2*2+1    'coin1, coun2 is unbiased coins. Use them as bits
   loop while sample>3             're-roll if got out of 1,2,3

   a(sample )=a(sample )+1
next
print a(1), a(2), a(3)


print "using Microsoft Linear congruential generator"
'from https://rosettacode.org/wiki/Linear_congruential_generator#Liberty_BASIC
global MSState

redim a(3)  'clear
for i= 1 TO N
SCAN
   sample=randMS() mod 3 +1
   a(sample )=a(sample )+1
next
print a(1), a(2), a(3)

'------------------------------------------------
function randMS()
' rand_n is in range 0 to 32767.
   MSState = (214013 * MSState + 2531011) mod (2 ^ 31)
   randMS = int(MSState / 2 ^ 16)
end function

function unBiasedCoin()
'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.
'Based on https://rosettacode.org/wiki/Unbias_a_random_generator
[again]
   n1 = rnd(1)>0.5
   n2 = rnd(1)>0.5
   if n1=n2 then [again]
   unBiasedCoin = n1
end function

EDIT
I believe line
sample = (temp mod 3) + 1
is wrong.

[tenochtitlanuk]

Confirm what Anatoly and B+ have said.
I've looked intensively at the small bias, including chi squared analysis and looking at frequency of digits and digit pairs.
For use in games LB's PRNG is more than good enough.

I use the re-seeding 'randomize' command to make EOR code/ciphering programs.

For use in mathematical number analysis it is less than perfect- for example using it to calculate pi will give a systematic slightly-wrong result.

The best ways are to-
-use the Mersenne generator, neatly packaged as a dll by Chris Iverson OR
-download a block of random digits from one of the free websites OR
-(expensive) buy a hardware random number generator ( & trust the manufacturer) OR
-(complicated unless you are an electronics whizz) follow instructions to make a RNG that uses something like 'noise' on a WiFi dongle.

[tsh73]

-use the Mersenne generator, neatly packaged as a dll by Chris Iverson

This one would work only for LibertyBASIC - bigger (paid) version of JB with DLLs and API support.

[B+]

Here is bplus random generator (using JB's) for numbers 1, 2, 3:

I am pleasantly surprised to get different distributions with this "3 card shuffle" for each run:

'r123 test.txt for JB 2 by bplus 2018-05-08

'initailize array
for i = 1 to 3
   r(i) = i
next

'check out 100
for i = 1 to 100
   print r123();" ";
   if i mod 20 = 0 then print
next

'test distribution
for i= 1 to 100000
   r = r123()
   d(r) = d(r) + 1
next
for i = 1 to 3
   print d(i), 100*d(i)/100000;"%"
next
end

function r123() 'shuffle and draw
   for i = 3 to 2 step -1
       r = int(rnd(0) * 3) + 1
       t = r(i)
       r(i) = r(r)
       r(r) = t
   next
   r123 = r(2) 'or r(1) or r(3)
end function


[tsh73]

What exactly surprises you?
They are random, so it'll never be exactly 33.33%
*waiting for 10 000 000 to finish*
Looks really close to me.
3330678       33.30678%
3333131       33.33131%
3336191       33.36191%

EDIT somehow initial post changes while I waited.
Ah well. Works nice anyway

[B+]

May 8, 2018 11:50:51 GMT tsh73 said:

What exactly surprises you?
They are random, so it'll never be exactly 33.33%
*waiting for 10 000 000 to finish*
Looks really close to me.
3330678       33.30678%
3333131       33.33131%
3336191       33.36191%

EDIT somehow initial post changes while I waited.
Ah well. Works nice anyway

I am surprised because I keep forgetting of rnd(0), the (0) part is NOT a seed.

I am not surprised it works because the shuffle was taught to me by you!

[bluatigro]

this can be useful

''bluatigro 22 apr 2018
''rnd module

function rnd.range( low , high )
 rnd.range = rnd(0) * ( high - low ) + low
end function
function rnd.int( low , high )
 rnd.int = int( rnd.range( low , high + 1 ) )
end function
function rnd.normal( mean , d )
 rnd.normal = sqr( abs( 2 * log( rnd(0) ) ) _
 * cos( rnd(0) * rnd.pi * 2 ) * d + mean
end function
function rnd.exp( mean )
 rnd.exp = log( rnd(0) ) * mean
end function



[bluatigro]

i tested the rnd function whit the following code

dim p( 4*4*4 )
global rock , paper , sisors , unkonwn , m1 , m2 , m3
global win , lose , tie
unknown = 0
rock = 1
paper = 2
sisors = 3
tel = 0
while tel < 3000
 c = AImove()
 h = rnd.move()
''  print "AI : " ; word$( "rock paper sisors" , c )
 if h = c then
   tie = tie + 1
''    print "tie ."
 else
   select case h
     case rock
       if c = paper then
         lose = lose + 1
''          print "You lose !"
       else
         win = win + 1
''          print "You win !!"
       end if
     case paper
       if c = sisors then
         lose = lose + 1
''          print "You lose !"
       else
         win = win + 1
''          print "You win !!"
       end if
     case else ''sisors
       if c = rock then
         lose = lose + 1
''          print "You lose !"
       else
         win = win + 1
''          print "You win !!"
       end if
   end select
 end if
''  input "again ? [ y / n ] : " ; yn$
 print tel
 tel = tel + 1
wend
print "Game over :"
print "played = " ; win + lose + tie
print "win = " ; win ; "   lose = " ; lose
print "points = " ; win - lose
end
function in( a , b , c )
 in = a + b * 4 + c * 16
end function
function humanmove()
 move$ = "5"
 while move$ < "1" or move$ > "3"
   print "1 = rock . 2 = paper . 3 = sisors ."
   input "your move " ; move$
 wend
 p( in( m1 , m2 , m3 ) ) = val( move$ )
 m1 = m2
 m2 = m3
 m3 = val( move$ )
 humanmove = m3
end function
function rnd.move()
 rnd.move = int( rnd(0) * 3 ) + 1
end function
function AImove()
 select case p( in( m1 , m2 , m3 ) )
   case rock
     uit = paper
   case paper
     uit = sisors
   case sisors
     uit = rock
   case else
     uit = int( rnd(0) * 3 ) + 1
 end select
 AImove = uit
end function

conclusion :
the ai can not find a patern

[zzz000abc]

one way to genreate random numbers write a function by taking k numbers.
we have to give equal chance to every numer and all k numbers shd not follow the same sequence.
so what we do we make all permutations of the given numbers.
once all permutations exhausted same permutations will repeat. in this way our randome function generates baised numbers.
other way is consider some irrational numbers and select one basing on some seed. now the digits of the iirational number becomes random numbers though every time same kind of baised ness is not observed.

[casco]

May 7, 2018 14:57:10 GMT combike said:

 
My solution: I decided to use the first digit after the decimal point of the random number MOD 3 to get the values I needed, and because the MOD operation will result in an unbalance of results, I eliminate all the 9s.
The code:
    do
      initial = rnd(1)
      temp = initial * 10
      sample = int(temp)
      sample = (temp mod 3) + 1
   loop while initial >= 0.9
This results in a more even distribution.
The question: Does anyone see anything wrong with what I've done, or is the RND() function kind of skewed?

Unfortunately your adjustment doesn't remedy to the problem, which exists with most of the pseudorandom number generators when they are subjected to sequences above a half of a million random numbers long. The local JB rnd() favors low random numbers and disfavors the high ones in the 0 to 1 range, and the bias behaves in almost linear manner. That means, if you generate integers in the 1-to-3 range, the sequence {1,2,3} would have the highest incidence, and the sequence (3,2,1} the lowest. I ran 60 samples, each of them containing 1,000,000 random numbers in the 1 - 3 range using your adjustment, and this is the result:

123______34

213______15

132_______9

231_______1

312_______1

321_______0


The chi-squared test soundly rejected the idea that the above distribution is due to the chance alone.

In the case of very extensive random numbers generating, I would recommend using the decimal digits of Pi as a source of random numbers. Being an irrational number, Pi is not made of self-repeating sequences and so it's bias free.

 Here is the code I used to check your rnd() adjustment. You can use it to test other pseudorandom number generators:


print "Progress..."



for samp=1 to 60 'number of samples

    dim a(3)

    string$=""



for size=1 to 1000000 'sample size

  r='RNG     <<<-----RNG comes here------

  a(r)=a(r)+1

next size



  b(1)=a(1)

  b(2)=a(2)

  b(3)=a(3)



   for i=1 to 2

   for j=i+1 to 3

     if a(i)<a(j) then

       hold=a(i)

       a(i)=a(j)

       a(j)=hold

     end if

   next j

   next i



   for i=1 to 3

   for j=1 to 3

     if a(i)=b(j) then string$=string$+str$(j)

   next j

   next i



   select case string$

    case "123": t123=t123+1

    case "132": t132=t132+1

    case "213": t213=t213+1

    case "231": t231=t231+1

    case "312": t312=t312+1

    case "321": t321=t321+1

   end select



    cls:print "sample ";samp;" done!"



next samp



cls

print 123, t123

print 132, t132

print 213, t213

print 231, t231

print 312, t312

print 321, t321



notice "end"

end