brianbec's WebLog

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 3:38 PM

so you are saying that you think a 4% error is large for a very basic PRNG?


seems normal to me.... how many values before you see a pattern (sequence repeats)?

I have to admit I have not been studying this recently but I have always understood that any "Random number" function is fake (Pseudo)
and when very simple tend to be very weak.

for example:
how are you seeding "Z" ?
and what is the period at which the function will start repeating values from a fixed seed.

as a simple example ( and forgive me if this is too basic etc...)

say the period was 1,000,000 values with a fixed seed number.
say you generated 3,500,000 values in each of say 20 runs.
if the not-so-random numbers tended to generate 1% more zeros in a run of one million then each run of 3.5 mill will show that 1% 2.5 times more than the 1's
but if the period was say a billion then the runs would seem much more random.

Denny Figuerres

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 4:01 PM

Sounds to me like your coin is damaged, or your tossing is biased.

If you convert the 31 bit number to a floating-point value between 0 and 1, and then check if it's less than 0.5, then most likely there's something wrong with the conversion.

If you test whether the current seed value falls in the first or second half of the range, then most likely your midpoint is off. It should be 2^30-1.

If you test straight away if the current seed value is odd, you should get the expected roughly 50/50 split.

Jeffrey Sax

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 5:13 PM

The period is 2^31-2, which is about 2 billion. I seeded the algorithm with System.Random, which, I believe, reads the current time. I counted bits in the output (0 has 0 1-bits, 1 has 1 1-bit, 2 has 1 1-bit, 3 has 2, and so on); that's a different test from converting the output to float or from flipping the coin based on even and odd. The reason I counted bits is that there are other PRNGs, most notably the Mersenne Twister (which has a period of something like 2^19,000), which has a bit ratio of 0.5000 to at least five decimal places over a sample of the same size as mine (600-odd billion). Another reason to like bit counting is that a single call of the PRNG can yield 32 coin flips if the bits are random.

So, what I'm really looking for is independent corroboration of the bit count ratio in the "minimal standar" PRNG, (z*16807)mod(2^31-1).

brianbec

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 5:30 PM

correction -- my samples were of size 600 million, not 600 billion :). I could easily run a sample the size of the whole period, and, if the bit-count ratio holds up, it would mean that in all the integers between 1 and 2^31-2, inclusive both ends, there are 48 percent 1 bits and 52 percent 0 bits. Adding 1 more bit to it and sampling via an algorithm with a vastly longer period via Mersenne (whose outputs are in [1..2^32-1] yield 50-50 ratio? I'll look into that first thing in the morning.

brianbec

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 5:45 PM

Are you taking into account the fact that there are 31 bits per sample? My own test gives a ratio of 0.5 to 6 decimal places. My bit counting function looks as follows. I use that n (the sample) can't be 0, so there is always at least one 1-bit. (not sure how this will look formatted)

Int32 newHeads = 0;
do {
n = n & (n-1);
newHeads++;
} while (n != 0);
heads += newHeads;
tails += 31 - newHeads;

Jeffrey Sax

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 6:09 PM

BTW, the statement "in all the integers between 1 and 2^31-2, inclusive both ends, there are 48 percent 1 bits and 52 percent 0 bits" is correct if you use 32 bits to represent the integers. The exact ratio is 15.5 to 16.5 = 48.4375% to 51.5625%. The first bit is always 0.

If you use 64 bit integers, the ratio is 24% to 76%.

Jeffrey Sax

# re: An Invitation to Flip Me Off@ Sunday, August 29, 2004 11:54 PM

AHA! Jeffrey, I think you've found it. Thanks. I will double check in the morning.

brianbec

# re: An Invitation to Flip Me Off@ Monday, August 30, 2004 4:15 PM

Okie doke, I've verified Jeffrey's result: counting bits out of 31 rather than out of 32 gives the expected statistical 50-pct ratio -- BUT, another mini-mystery has surfaced. Let z be the sample, p=2^31-1, a=16807, and m=2^32-1. I would expect that

double d = ((a*z)mod p)/(p-1)
would be uniform (0.0 .. 1.0], so

d * m would have a statistical 50-pct ratio out of 32 bits, but it doesn't. Out of 100000 samples, it generates 51.0729. Consider the following Mathematica code... Explanations anyone?

The modulus:

In[1]:=
m=2^31-1

Out[1]=
2147483647

The multiplier:

In[2]:=
a=16807;

Verify that the multiplier will generate the entire group (see blogs for rationale):

In[3]:=
(m-1)==MultiplicativeOrder[a,m]

Out[3]=
True

The Pseudo-Random Number Generator (PRNG):

In[4]:=
p[z_]:=Mod[(z*a),m]

Unit test for the first few values:

In[5]:=
Module[{z=92},Table[z=p[z],{16}]]

Out[5]=
{1546244,217919144,1107435073,420503362,31322857,309764084,704599460,\
978294662,1063582802,2129759233,606000835,1688579771,963816092,387908923,\
1972400216,1572855220}

The bit-counting function. Takes three arguments: a function to iterate, a seed, and a number of samples to take. It counts the number of on bits (nh for "number of heads") and the number of off bits (nt for "number of tails"), but ASSUMES 32-bit words, as verified by the output. We shall expect approximately 0.484375 = 15.5/32 for nh on PRNG "p", because "p" produces 31-bit integers:

In[6]:=
bc[f_,s_,n_]:=Module[{z=s,nh=0,nt=0,nb=0},
Do[
z=f[z];
nh+=DigitCount[z,2,1];
nt+= DigitCount[z,2,0],{n}];
nb=nh+nt;
{nh,nt,N[nh*100/nb],N[nt*100/nb],N[nh/nt]}]

Test it on the PRNG, "p":

In[7]:=
bc[p,92,100000]

Out[7]=
{1550834,1449550,51.6879,48.3121,1.06987}

In[8]:=
p1[z_]:=IntegerPart[(N[(2^32-1)*p[z]/(m-1)])]

In[9]:=
bc[p1,92,100000]

Out[9]=
{1583561,1517029,51.0729,48.9271,1.04386}

brianbec

# re: An Invitation to Flip Me Off@ Monday, August 30, 2004 4:41 PM

The ratio of floating-point samples above 0.5 to those below shows acceptable statistics:

In[28]:=
uc[f_,s_,n_,divisor_]:=Module[{z=s,nhi=0,nlo=0,d=0.0},
Do[
z=f[z];
d=N[z/divisor];
If[d>0.5,nhi++,nlo++],{n}];
{nhi,nlo,N[nhi*100/n],N[nlo*100/n],N[nhi/nlo]}]

In[32]:=
uc[p,92,1000000,(m-1)]

Out[32]=
{500655,499345,50.0655,49.9345,1.00262}

brianbec