[FIXED] Drawing random numbers and advancing the bit generator yields different results

Issue

I am trying to advance my NumPy random number generator to mimic drawing n random numbers. However, when I compare the state of the bit generator after advancing n numbers with simply drawing n numbers, the results differ after some point.

I know that some random numbers might require to advance the state multiple times but I find it odd that the results match for the first 40-50 drawn numbers and then diverge.

I am trying to use advance because at some times I have to draw random numbers to keep reproducing the same results although for some use-cases I already know that I don't use the numbers. In this case, I want to simply (and fast) advance the random number generator but ensuring that all the numbers drawn remain the same.

How can I ensure reproducability and still use features like advance?

If I use a random number generator from numpy, I compare the state of the bit_generator for the cases when I advance 10 deltas at a time with drawing the same number of random variable (standard_normal) at a time. After few iterations, the rng.bit_generator.state diverges. I don't know why because I would expect it do diverge instantly or never.

rng = np.random.default_rng(2024)
s = rng.bit_generator.state

step = 10
num = 5

# Draw Random numbers using standard normal variable

rng.bit_generator.state = s

x = list()

for i in range(num):
    rng.standard_normal(step)
    x.append(rng.bit_generator.state['state']['state'])
    print('{i:3d} - {state:40d} {inc:40d}'.format(i=i, **rng.bit_generator.state['state']))

# Reset state and advance random number generator
rng.bit_generator.state = s

y = list()

for i in range(num):
    rng.bit_generator.state = rng.bit_generator.advance(step).state
    y.append(rng.bit_generator.state['state']['state'])
    print('{i:3d} - {state:40d} {inc:40d}'.format(i=i, **rng.bit_generator.state['state']))

The output looks like this for me. The results differ for iteration 2, i.e. the state is different.

  0 -  317319928805135160732659717497650841180  263843294879837360010514471918415607657
  1 -  152261761237047187080279356346465986730  263843294879837360010514471918415607657
  2 -   12941317612297752247903241541248335614  263843294879837360010514471918415607657
  3 -  167375039976652124381293896870728401724  263843294879837360010514471918415607657
  4 -   78535979019656695810700673044059278986  263843294879837360010514471918415607657

  0 -  317319928805135160732659717497650841180  263843294879837360010514471918415607657
  1 -  152261761237047187080279356346465986730  263843294879837360010514471918415607657
  2 -  262244773524709030991986015855808125320  263843294879837360010514471918415607657
  3 -  174718512225906865463277455530823048310  263843294879837360010514471918415607657
  4 -  213819718814773627130094102456555149556  263843294879837360010514471918415607657

Solution

Why are you assuming that advance does the same as generating a sample of a given type>

In [715]: rng=np.random.default_rng()

In [716]: timeit rng.bit_generator.advance(1000)
2.66 µs ± 5.72 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)

Uniform and normal take different times, suggesting something different in how they are generated:

In [717]: timeit rng.random(1000)
10.2 µs ± 9.81 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)

In [718]: timeit rng.standard_normal(1000)
19.3 µs ± 41.7 ns per loop (mean ± std. dev. of 7 runs, 10,000 loops each)

advance docs says:

Advancing a RNG updates the underlying RNG state as-if a given
number of calls to the underlying RNG have been made. In general
there is not a one-to-one relationship between the number output
random values from a particular distribution and the number of
draws from the core RNG.  

uniform is the simpler distribution, taking less time, and more likely to have a one-to-one relation with the advance:

In [727]: rng = np.random.default_rng(2024)
     ...: s = rng.bit_generator.state
    
In [728]: rng.bit_generator.state = s
     ...: for i in range(5):
     ...:     rng.bit_generator.state = rng.bit_generator.advance(10).state
     ...:     print('{i:3d} - {state:40d} {inc:40d}'.format(i=i, **rng.bit_generator.state['state']))
     ...:     
  0 -  317319928805135160732659717497650841180  263843294879837360010514471918415607657
  1 -  152261761237047187080279356346465986730  263843294879837360010514471918415607657
  2 -  262244773524709030991986015855808125320  263843294879837360010514471918415607657
  3 -  174718512225906865463277455530823048310  263843294879837360010514471918415607657
  4 -  213819718814773627130094102456555149556  263843294879837360010514471918415607657

In [730]: rng.bit_generator.state = s
     ...: for i in range(5):
     ...:     rng.standard_normal(10)
     ...:     print('{i:3d} - {state:40d} {inc:40d}'.format(i=i, **rng.bit_generator.state['state']))
     ...:     
  0 -  317319928805135160732659717497650841180  263843294879837360010514471918415607657
  1 -  152261761237047187080279356346465986730  263843294879837360010514471918415607657
  2 -   12941317612297752247903241541248335614  263843294879837360010514471918415607657
  3 -  167375039976652124381293896870728401724  263843294879837360010514471918415607657
  4 -   78535979019656695810700673044059278986  263843294879837360010514471918415607657

In [731]: rng.bit_generator.state = s
     ...: for i in range(5):
     ...:     rng.random(10)
     ...:     print('{i:3d} - {state:40d} {inc:40d}'.format(i=i, **rng.bit_generator.state['state']))
     ...:     
  0 -  317319928805135160732659717497650841180  263843294879837360010514471918415607657
  1 -  152261761237047187080279356346465986730  263843294879837360010514471918415607657
  2 -  262244773524709030991986015855808125320  263843294879837360010514471918415607657
  3 -  174718512225906865463277455530823048310  263843294879837360010514471918415607657
  4 -  213819718814773627130094102456555149556  263843294879837360010514471918415607657

ramdom does match the advance, standard_normal does not.

Answered By - hpaulj