Should mvnrnd Always Advance the State of the Global Stream
2 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Consider the following:
>> mu=[1 1]; Sigma=eye(2);
rng('default')
preu1 = rand(1,3);
n1 = mvnrnd(mu,Sigma);
u1 = rand(1,3);
rng('default')
preu2 = rand(1,3);
n2 = mvnrnd(mu,Sigma);
u2 = rand(1,3);
rng('default')
n3 = mvnrnd(mu,0*Sigma);
u3 = rand(1,3);
>> [isequal(u1,u2) isequal(u2,u3) isequal(u3,preu2)]
ans =
1×3 logical array
1 0 1
Apparently mvnrnd doesn't actually call randn if it detects that the input covariance is zero. That may be good for efficiency, but is it the best behavior for repeatability? This behavior seems to be a contradiction to the general direction of how to manage the Global Stream to reproduce results. doc mvnrnd is silent on how it handles this special case.
0 Commenti
Risposte (1)
Steven Lord
il 26 Mag 2020
There's no requirement that these two lines of code that call the same function with different inputs need follow identical code paths.
myfun(A, B)
myfun(A, C)
The lines of code in your example look similar to the human eye, but I could rewrite them to be less similar but equivalent to what you wrote.
x1 = mvnrnd(mu, eye(2))
x2 = mvnrnd(mu, zeros(2))
9 Commenti
Steven Lord
il 6 Lug 2023
randn isn't free. But it's a price I'd be willing to pay to make it easier to have repeatable Monte Carlo simulations.
Is every user of mvnrnd also willing to pay that price? I'm fairly sure the answer to that question is no.
Scanning back through the three year history of this discussion, I noted this section that I don't think I paid close enough attention to back when it was posted.
As it stands, consider an application with a sequence of calls to various RNGs, including mvnrnd, and assume that anything that controls the dimensions of the outputs of the RNGs is fixed. However, the values of the input parameters that are used to call the RNG functions are provided by the user. Would it not be reasonable to expect that executing the sequence twice, and executing rng('default') prior each execution, result in the same outputs for those RNGs that used the same parameters in both runs? As shown in the example, if on the second run the user sets Sigma to 0 for a call to mvnrnd, then all rand* calls downcode from that call to mvnrnd that use the Global Stream will be affected as well.
Your expectation is reasonable. The key point is that for the same user inputs the code should give the same results. But in your example, Sigma is one of those user inputs! So in the following example I would expect x0 and x1 to be the same, but I would not necessarily expect x2 to be the same as x0 or x1 because they're using different values of Sigma. Since mySimulation calls mvnrnd internally then calls some other random number generator function (let's say randi), the call that generates x0 may get different results internally from randi than the call that generates x2.
Sigma = 0;
rng default
x0 = mySimulation(Sigma);
rng default
x1 = mySimulation(Sigma);
Sigma = 1;
rng default
x2 = mySimulation(Sigma);
isequal(x0, x1)
isequal(x0, x2)
I'm pretty sure that I've even seen somehwere in the doc that advises to not do things like this, because it shouldn't be needed.
I'm pretty sure you're referring to the Note on this documentation page. That's intended to warn against doing something like the following (commented out since calling rng shuffle 1e5 times takes long enough that the MATLAB session MATLAB Answers uses to run code reaches the time limit.) Shuffling the random number generator 1e5 times doesn't necessarily make the numbers "more random" like the comment states the code author believes.
But reseeding the generator once before each time you run your simulation (for reproducibility) as I did above in the x0 / x1 / x2 example is fine.
%{
n = 1e5;
x1 = zeros(1, n);
for k = 1:n
rng shuffle % To make the numbers "more random"
x1(k) = rand;
end
%}
function y = mySimulation(Sigma)
% "Burn off" Sigma + 10 numbers then return a random integer in the range [1, 10]
randn(1, Sigma+10);
y = randi(10, 1);
end
Vedere anche
Categorie
Scopri di più su MATLAB Mobile Fundamentals in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!