Menu
More specifically, numpy:
In [24]: a=np.random.RandomState(4)
In [25]: a.rand()
Out[25]: 0.9670298390136767
In [26]: a.get_state()
Out[26]:
('MT19937',
array([1248735455, ..., 1532921051], dtype=uint32),
2,0,0.0)
octave:
octave:17> rand('state',4)
octave:18> rand()
ans = 0.23605
octave:19> rand('seed',4)
octave:20> rand()
ans = 0.12852
Octave claims to perform the same algorithm (Mersenne Twister with a period of 2^{19937-1})
Anybody know why the difference?
212
The random Module Probably the most widely known tool for generating random data in Python is its random module, which uses the Mersenne Twister PRNG algorithm as its core generator. If you run this code yourself, I'll bet my life savings that the numbers returned on your machine will be different.
NumPy random for generating an array of random numbers 10 000 calls, and even though each call takes longer, you obtain a numpy. ndarray of 1000 random numbers. The reason why NumPy is fast when used right is that its arrays are extremely efficient. They are like C arrays instead of Python lists.
Python defines a set of functions that are used to generate or manipulate random numbers through the random module. Functions in the random module rely on a pseudo-random number generator function random(), which generates a random float number between 0.0 and 1.0.
Indeed, whenever we call a python function, such as np. random. rand() the output can only be deterministic and cannot be truly random. Hence, numpy has to come up with a trick to generate sequences of numbers that look like random and behave as if they came from a purely random source, and this is what PRNG are.
Unfortunately the MT19937 generator in Octave does not allow you to initialise it using a single 32-bit integer as np.random.RandomState(4) does. If you use rand("seed",4) this actually switches to an earlier version of the PRNG used previously in Octave, which PRNG is not MT19937 at all, but rather the Fortran RANDLIB.
It is possible to get the same numbers in NumPy and Octave, but you have to hack around the random seed generation algorithm in Octave and write your own function to construct the state vector out of the initial 32-bit integer seed. I am not an Octave guru, but with several Internet searches on bit manipulation functions and integer classes in Octave/Matlab I was able to write the following crude script to implement the seeding:
function state = mtstate(seed)
state = uint32(zeros(625,1));
state(1) = uint32(seed);
for i=1:623,
tmp = uint64(1812433253)*uint64(bitxor(state(i),bitshift(state(i),-30)))+i;
state(i+1) = uint32(bitand(tmp,uint64(intmax('uint32'))));
end
state(625) = 1;
Use it like this:
octave:9> rand('state',mtstate(4));
octave:10> rand(1,5)
ans =
0.96703 0.54723 0.97268 0.71482 0.69773
Just for comparison with NumPy:
>>> a = numpy.random.RandomState(4)
>>> a.rand(5)
array([ 0.96702984, 0.54723225, 0.97268436, 0.71481599, 0.69772882])
The numbers (or at least the first five of them) match.
Note that the default random number generator in Python, provided by the random module, is also MT19937, but it uses a different seeding algorithm, so random.seed(4) produces a completely different state vector and hence the PRN sequence is then different.
72
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With
