Abhinav Gupta
Abhinav Gupta

Reputation: 207

matlab random number generation in parfor loop

I want to generate same normal random numbers for loop and parfor loop. Following MATLAB documentation I tried three different method:

Method 1: Using rng(seed, 'twister')

N = 1;

for ind = 1:10
    rng(ind, 'twister')
    samps(ind) = normrnd(0,1,N,1);
end

plot(samps); hold on

parfor ind = 1:10
    rng(ind, 'twister')
    samps(ind) = normrnd(0,1,N,1);
end

scatter(1:length(samps), samps)
legend({'using for loop', 'using parfor loop'})

By using rng twister method

Method 2: Using RandStream Method

N = 1;
 
sc = parallel.pool.Constant(RandStream('Threefry'));
for ind = 1:10
    stream = sc.Value;
    stream.Substream = ind;
    samps(ind) = normrnd(0,1,N,1);
end

plot(samps); hold on

sc = parallel.pool.Constant(RandStream('Threefry'));
parfor ind = 1:10
       stream = sc.Value;
       stream.Substream = ind;
       samps(ind) = normrnd(0,1,N,1);
end

scatter(1:length(samps), samps)
legend({'using for loop', 'using parfor loop'})

Using Threefry stream-substream method

Method 3: Using another RandStream method

N = 1;
 
stream = RandStream('mrg32k3a');
for ind = 1:10
    stream.Substream = ind;
    samps(ind) = normrnd(0,1,N,1);
end

plot(samps); hold on

stream = RandStream('mrg32k3a');
parfor ind = 1:10
       set(stream,'Substream',ind);
       samps(ind) = normrnd(0,1,N,1);
end

scatter(1:length(samps), samps)
legend({'using for loop', 'using parfor loop'})

Another randStream method

My question is why the last two methods are not working. If I understood MATLAB documentation correctly, they should have worked. Also, is there anything wrong with using rng(seed, 'twister') method with different seeds to produce statistically independent samples?

Edit: Here is the link to the MATLAB documentation that I used.

Upvotes: 1

Views: 679

Answers (3)

Edric
Edric

Reputation: 25160

As @Cris points out, normrnd uses the current global stream. Here's how to make your 3rd case match between client and workers - you need to set the global stream explicitly. (Note that I'm reverting the global stream immediately after using it too).

N = 1;

stream = RandStream('mrg32k3a');
for ind = 1:10
    stream.Substream = ind;
    prev = RandStream.setGlobalStream(stream);
    samps(ind) = normrnd(0,1,N,1);
    RandStream.setGlobalStream(prev);
end

stream = RandStream('mrg32k3a');
parfor ind = 1:10
    set(stream,'Substream',ind);
    prev = RandStream.setGlobalStream(stream);
    pf_samps(ind) = normrnd(0,1,N,1);
    RandStream.setGlobalStream(prev);
end

With this code, samps and pf_samps are identical.

Also, you should not use "twister" in parallel, since the random samples you obtain by setting the seed in this way will not have as good qualities as if you use a "parallel" generator like mrg32k3a. (I'm no expert, but my understanding is that the twister "seed" values are not able to select statistically-independent sequences of samples in the way that streams/substreams can)

In recent versions of MATLAB, faster generators that support multiple streams are available. See the doc https://www.mathworks.com/help/matlab/math/creating-and-controlling-a-random-number-stream.html#brvku_2 . Philox and Threefry in particular are good parallel generators.

Upvotes: 2

Cris Luengo
Cris Luengo

Reputation: 60761

When you do normrnd(0,1,N,1), you are not using stream to generate a random number, but rather the default RNG. Instead, use randn(stream,N,1), see the documentation.

As far as I can tell, normrnd and other generators from the Statistics and a Machine Learning toolbox always use the default RNG, and cannot use a custom RandStream object.

Upvotes: 0

Nigel
Nigel

Reputation: 3291

I'm inferring that your concern is that the pseudorandomly generated numbers for your for and parfor loops did not match. This is however expected behavior. You used a MRG32K3A pseudorandom algorithm, which "does not make a guarantee that a fixed seed and a fixed series of calls to mrg32k3a.RandomState methods using the same parameters will always produce the same results." Source.

So it appears that everything is working correctly.

Upvotes: 1

Related Questions