Randperm is not random! How to REALLY randomise each randperm() call?
in [Matlab]
|
From: David Painter on 23 Feb 2010 06:46 Hi, I've finally realised randperm is subject to the same "limitations" as the rand functions as it contains rand (obviously) and so operates on the default random stream. :'-( This may seem obvious, but it should be documented explicitly in the help file. No warning/hint is given in Matlab R2008B. Wondering why it is not advisable to seed the random by the block more than once per Matlab session? s = RandStream.create('mt19937ar','seed',sum(100*clock)); RandStream.setDefaultStream(s); randperm(10) s = RandStream.create('mt19937ar','seed',sum(100*clock)); RandStream.setDefaultStream(s); randperm(10) I want my randperms to produce different results on each run. How can I do it? Thanks, David
From: John D'Errico on 23 Feb 2010 06:55 "David Painter" <david.ross.painter(a)gmail.com> wrote in message <hm0f5s$ct6$1(a)fred.mathworks.com>... > Hi, > > I've finally realised randperm is subject to the same "limitations" as the rand functions as it contains rand (obviously) and so operates on the default random stream. :'-( This may seem obvious, but it should be documented explicitly in the help file. No warning/hint is given in Matlab R2008B. > > Wondering why it is not advisable to seed the random by the block more than once per Matlab session? > > s = RandStream.create('mt19937ar','seed',sum(100*clock)); > RandStream.setDefaultStream(s); > randperm(10) > s = RandStream.create('mt19937ar','seed',sum(100*clock)); > RandStream.setDefaultStream(s); > randperm(10) > > I want my randperms to produce different results on each run. How can I do it? Do NOT reseed the random generator multiple times within one session. This will make the random number generator LESS random, not more so. And randperm DOES return different results on each call. randperm(10) ans = 6 3 7 8 5 1 2 4 9 10 randperm(10) ans = 6 1 7 4 9 5 8 3 10 2 randperm(10) ans = 2 10 8 9 1 5 7 6 3 4 John
From: Jan Simon on 23 Feb 2010 10:30 Dear David! > s = RandStream.create('mt19937ar','seed',sum(100*clock)); > RandStream.setDefaultStream(s); > randperm(10) > s = RandStream.create('mt19937ar','seed',sum(100*clock)); > RandStream.setDefaultStream(s); > randperm(10) If you believe that mt19937ar is a valid random number generator, there is no need to seed it twice. If you don't believe that mt19937ar produces valid random numbers, use another generator. I personally do trust Matlab's random number generator and I'm convinced that RANDPERM is not biased. As far as I understood the only drawback of RANDPERM is the stable SORT: If a random value appears twice, the order of these elements is kept. Fortunately the chance to get two equal values in a RAND vector can be neglected. Another drawback of RANDPERM is the speed for large arrays, which is limited by the sorting. The Fisher-Yates shuffle is simply faster (60% of RANPERM time for 10 elements, 46% of RANDPERM for 1E6 elements): function X = RandInplacePermute(X) for i = 1:numel(X) w = ceil(rand * i); t = X(w); X(w) = X(i); X(i) = t; end The main difference appears for large arrays, because RANDPERM needs 4 temporary vectors of the same size (see source of RANDPERM): 1. RAND(1, n) 2. SORT creates one data copy for the ignored sorted output and one for the indices. 3. X = X(randperm(numel(X)) creates a further copy. I think this could help to re-use the memory: X(:) = X(randperm(numel(X)) Anyhow, I think RandInplacePermute is too trivial to be published on the FEX - different opinions? Kind regards, Jan
From: David Painter on 23 Feb 2010 18:19 "John D'Errico" <woodchips(a)rochester.rr.com> wrote in message <hm0fmq$ft2$1(a)fred.mathworks.com>... > "David Painter" <david.ross.painter(a)gmail.com> wrote in message <hm0f5s$ct6$1(a)fred.mathworks.com>... > > Hi, > > > > I've finally realised randperm is subject to the same "limitations" as the rand functions as it contains rand (obviously) and so operates on the default random stream. :'-( This may seem obvious, but it should be documented explicitly in the help file. No warning/hint is given in Matlab R2008B. > > > > Wondering why it is not advisable to seed the random by the block more than once per Matlab session? > > > > s = RandStream.create('mt19937ar','seed',sum(100*clock)); > > RandStream.setDefaultStream(s); > > randperm(10) > > s = RandStream.create('mt19937ar','seed',sum(100*clock)); > > RandStream.setDefaultStream(s); > > randperm(10) > > > > I want my randperms to produce different results on each run. How can I do it? > > Do NOT reseed the random generator multiple times > within one session. This will make the random > number generator LESS random, not more so. > > And randperm DOES return different results on each > call. > > randperm(10) > ans = > 6 3 7 8 5 1 2 4 9 10 > > randperm(10) > ans = > 6 1 7 4 9 5 8 3 10 2 > > randperm(10) > ans = > 2 10 8 9 1 5 7 6 3 4 > > John Thanks to everyone for their helpful replies. I am guilty of posting while tired! John, your're right. Randperm produces different results on each call within a session (open Matlab, exit Matlab). The issue is how to have different sequences for different sessions. One call to seed the stream based on the clock is sufficient because clock is likely to be different on each call to: s = RandStream.create('mt19937ar','seed',sum(100*clock)) Is that right?
From: Greg Heath on 25 Feb 2010 04:20 On Feb 24, 11:07 am, Peter Perkins <Peter.Perk...(a)MathRemoveThisWorks.com> wrote: > On 2/23/2010 6:19 PM, David Painter wrote: > > > One call to seed the stream > > based on the clock is sufficient because clock is likely to be different > > on each call to: > > > s = RandStream.create('mt19937ar','seed',sum(100*clock)) > > > Is that right? > > It is, although > > s = RandStream('mt19937ar','seed',sum(100*clock)) > > s marginally shorter to type (since you're creating only one random stream), and don't forget that you need to make your new stream the default random number stream, > > RandStream.setDefaultStream(s) > > You may find that simply reseeding the existing default stream > > reset(RandStream.getDefaultStream,sum(100*clock)) > > at the beginning of each MATLAB session (e.g., in your startup file) is even easier. For today (Feb 25, 2010) min(sum(clock)) = 2010 + 2 + 25 + 0 + 0 + 0 = 2037 max(sum(clock)) = 2010 + 2 + 25 + 23 + 59 + 59 = 2178 Therefore, Jseed = 100*sum(clock) only allows 100*(23+59+59) = 14,100 integer seeds per day whereas the total number of available integer seeds is 4,294,967,296 over the interval [ 0, 2^32-1] = [0 4,294,967,295] Perhaps Jseed = 3e7*sum(clock-clock0) where clock0 = 2010+2+25 (depends on the data) ( or one of many other variations) might yield a better diversity of random streams. Hope this helps. Greg
|
Pages: 1 Prev: true color compression image Next: Using csaps Cubic Smoothing spline |