What RNGkind does (was Re: Bug in rnorm. (PR#1664))
ripley@stats.ox.ac.uk
ripley@stats.ox.ac.uk
Fri, 14 Jun 2002 10:03:45 +0100 (BST)
Just so we understand:
> RNGkind()
[1] "Marsaglia-Multicarry" "Kinderman-Ramage"
> .Random.seed <- c(0, 1, 2, 3)
> RNGkind()
[1] "Marsaglia-Multicarry" "Kinderman-Ramage"
> runif(1)
[1] 0.03381877
> RNGkind()
[1] "Wichmann-Hill" "Kinderman-Ramage"
?RNGkind says
`RNGkind' returns a two-element character vector of the RNG and
normal kinds in use before the call,
^^^^^^
However, any use of the random number generators themselves loads
.Random.seed (if it exist), and that will set the types of random
numbers in use. So if a user changes .Random.seed (for example by
loading a saved workspace) the types will change. RNGkind() reports
retrospectively (as it is documented to), not prospectively.
Now we could alter RNGkind() to read in the current .Random.seed (if there
is one) although it is a bit tricky to decide what to do if the value is
invalid. (The random-number generation code throws an error and
expects thse user to fix, e.g remove, .Random.seed.)
On 14 Jun 2002, Peter Dalgaard BSA wrote:
> Duncan Murdoch <dmurdoch@pair.com> writes:
>
> > On Fri, 14 Jun 2002 00:39:58 +0200 (MET DST), p.dalgaard@biostat.ku.dk
> > wrote:
> >
> > >With your script, I get a nice horizontal line both with 1.4.0 and
> > >1.5.0-patched.
> > >
> > >(RedHat Linux 7.1 on a PC)
> >
> > I get what Rolf described in 1.5.0-patched for Windows with the
> > default generators:
> >
> > RNGkind()
> > [1] "Marsaglia-Multicarry" "Kinderman-Ramage"
> >
> > I switched to the Ahrens-Dieter normal RNG (using
> > RNGkind(,'Ahrens-Dieter')), and things were fine.
> >
> > Keeping the default normal generator but switching to the
> > Mersenne-Twister uniform also fixed it:
> >
> > > RNGkind()
> > [1] "Mersenne-Twister" "Kinderman-Ramage"
> >
> > Is it possible you're not using the default generators? If you are,
> > this is quite weird, because we're on the same hardware. Why would
> > the OS or compiler matter in a calculation like this??
> >
> > Duncan Murdoch
> >
>
> Now this is pretty darn odd....
>
> This is what happened to me yesterday with 1.4.0:
>
> [Previously saved workspace restored]
>
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 39
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 37
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 63
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 50
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 51
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 51
>
> This morning, I get
>
> > RNGkind()
> [1] "Marsaglia-Multicarry" "Kinderman-Ramage"
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 27
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 23
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 29
> > m <- sapply(1:1000, function(i)max(rnorm(1500)^2))
> > sum( 1 - pchisq(m,1)^1500 < .05)
> [1] 29
>
> with the SAME version of R. The only difference is the startup
> directory, and hence the saved workspace. This contains:
>
> > ls(all=T)
> [1] ".Random.seed" "x" "y"
> > .Random.seed
> [1] 0 25300 11635 10783
>
> ..whereas 1.5.0 which now again shows trouble has
>
> > .Random.seed
> [1] 1 196857153 319227924
>
> and 1.4.0 in a clean dir has this after a few iterations of the above
>
> > .Random.seed
> [1] 1 1297109091 47580530
>
>
> Oho! I think I get it: the "1" in the seed actually *defines* the RNG,
> no matter what RNGkind() is telling me??
--
Brian D. Ripley, ripley@stats.ox.ac.uk
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272860 (secr)
Oxford OX1 3TG, UK Fax: +44 1865 272595
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._