Is there any way to generate uncorrelated random variables using Python?

Question

Suppose I want to generate two random variables X and Y which are uncorrelated and uniformly distributed in [0,1].

The very naive code to generate such is the following, which calls the random function twice:

import random 
xT=0 
yT=0 
xyT=0 
for i in range(20000):
    x = random.random()
    y = random.random()
    xT += x
    yT += y
    xyT += x*y

xyT/20000-xT/20000*yT/20000

However, the random number is really a pseudo-random number which is generated by a formula, therefore they are correlated.

How to generate two uncorrelated (or correlated as little as possible) random variables?

Answer 1

The math on RNGs is solid. These days most popular implementations are too. As such, your conjecture of

is generated by a formula, therefore they are correlated.

is incorrect.

But if you really truly deeply think that way, there is an out: hardware random number generators. The site at random.org has been providing hardware RNG draws "as a service" for a long time. Here is an example (in R, which I use more, but there is an official Python client):

R> library(random)
R> randomNumbers(min=1, max=20000)    # your range, default number
         V1    V2    V3    V4    V5
 [1,]   532 19452  5203 13646  5462
 [2,]  4611 10814  3694 12731   566
 [3,] 11884 19897  1601 10652   791
 [4,] 17427  9524  7522  1051  9432
 [5,]  5426  5079  2232  2517  4883
 [6,] 13807  9194 19980  1706  9205
 [7,] 13043 16250 12827  2161 10789
 [8,]  7060  6008  9110  8388  1102
 [9,] 12042 19342  2001 17780  3100
[10,] 11690  4986  4389 14187 17191
[11,] 19574 13615  3129 17176  5590
[12,] 11104  5361  8000  5260   343
[13,]  7518  7484  7359 16840 12213
[14,] 14914  1991 19952 10127 14981
[15,] 13528 18602 10182  1075 16480
[16,]  9631 17160 19808 11662 10514
[17,]  4827 13960 17003   864 11159
[18,]  8939  7095 16102 19836 15490
[19,]  8321  6007  1787  6113 17948
[20,]  9751  7060  8355 19065 15180
R> 

Edit: The OP seems unconvinced, so there is a quick reproducible simulation (again, in R because that is what I use):

R> set.seed(42)               # set seed for RNG
R> mean(replicate(10, cor(runif(100), runif(100))))
[1] -0.0358398
R> mean(replicate(100, cor(runif(100), runif(100))))
[1] 0.0191165
R> mean(replicate(1000, cor(runif(100), runif(100))))
[1] -0.00117392
R> 

So you see that as we go from 10 to 100 to 1000 replications of just 100 U(0,1), the correlations estimate goes to zero.

We can make this a little nice with a plot, recovering the same data and some more:

R> set.seed(42)
R> x <- 10^(1:5)   # powers of ten from 1 to 5, driving 10^1 to 10^5 sims
R> y <- sapply(x, function(n) mean(replicate(n, cor(runif(100), runif(100)))))
R> y    # same first numbers as seed reset to same start
[1] -0.035839756  0.019116460 -0.001173916 -0.000588006 -0.000290494
R> plot(x, y, type='b', main="Illustration of convergence towards zero", log="x")
R> abline(h=0, col="grey", lty="dotted")

Answer 2

Short answer: Use a Bays-Durham shuffle on your random seeds.

Longer answer:

I'm sure you know that the pseudo-random numbers given by computer algorithms are not truly random--they are just meant to pass most tests for randomization and thus be "good enough" for most uses. The same is true for uncorrelated random variables: you will never get truly uncorrelated ones, but your goal should be to get them to pass as many correlation tests as possible and be "good enough" for your purposes.

The main way that the standard, linear congruential modulators fail correlation tests is when you look at a small region of the 2-space generated by the numbers. The pairs of numbers show obvious lines when graphed and thus are not truly uncorrelated. This only matters when you look at a very small region of all the generated number pairs. Is this what you need to do? Note that Python's random() function uses the "Mersenne Twister" rather than a linear congruential modulator, so it is less likely to fail such a test. See Wikipedia's list of the disadvantages of the Mersenne Twister to see if Python's random number generator is or is not good enough for your purposes. Note that Python's implementation is shown in detail later in the page.

I wrote routines in Borland Delphi (Object Pascal and x86 assembler) to avoid correlation. I have switched to Python but have not yet rewritten those routines. The idea of a Bays-Durham shuffle is to use the built-in random number generator to give you a random integer (the one that is used to make a floating-point number between 0 and 1). You then use that integer to point into an array of previously-generated random integers. You pick the integer pointed to and replace it in the array with your newly generated integer. Return the integer that used to be in the array, and if need be convert that to a number between 0 and 1.

I implemented this with an array of 32 integers and tested this new generator. This now passed the correlation test that Delphi's random number generator failed. I repeat, this would not pass all correlation tests, but it did pass more than the standard generator, and it was definitely good enough for my uses.

If you need to see a Python implementation of this, ask and I'll try to make the time to write one. Until then, look up "Bays-Durham shuffle". I learned of it from the book Numerical Recipes. Here is a Fortran version of that chapter. Here is the entire 2nd edition in C in Empanel format and here it is in PDF--find Chapter 7 Section 7.1. Source code for out-of-date editions is available in a variety of languages, including Fortran (I think), C, and Pascal. I downloaded the 2nd edition C version text and 1st edition Pascal code some years ago and used those for my own coding in Pascal.