Bug in rnorm. (PR#1664)
Peter Dalgaard BSA
p.dalgaard@biostat.ku.dk
14 Jun 2002 10:14:29 +0200
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??
--
O__ ---- Peter Dalgaard Blegdamsvej 3
c/ /'_ --- Dept. of Biostatistics 2200 Cph. N
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard@biostat.ku.dk) FAX: (+45) 35327907
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._