[R] random sequence

[Prof Brian D Ripley]

On Mon, 26 Apr 1999, Paul Gilbert wrote:

> Brian
> 
> >What exactly do you want?  My class exercise is to implement the
> >Wichmann-Hill algorithm in C/Fortran and use it with S-PLUS! But
> >while one is at it, why not use a decent generator with both?
> 
> It would certainly be nice to have all the R generators in S, but for now I
> would be happy to have just one. What I would like is to be able to run the
> programs which test my code, and get the same random sequences so that I don't
> need different check values for R and Splus.

>From my files (from before R was thought of):

WichHill <- function(n)
{
        output <- numeric(n)
        seed <- get(".WHseed", frame=0)
        x <- seed[1]; y <- seed[2]; z <- seed[3]
        for(i in 1:n){
                x <- (171*x) %% 30269
                y <- (172*y) %% 30307
                z <- (170*z) %% 30323
                output[i] <- (x/30269 + y/30307 + z/30323) %% 1.0
        }
        assign(".WHseed", c(x,y,z), frame=0, immediate=T)
        output
}

> assign(".WHseed",  1:3, frame=0)
> WichHill(10)
 [1] 0.03381877 0.77754189 0.05273525 0.74462407 0.49036219 0.98285437
 [7] 0.80915099 0.71338138 0.80102091 0.98958603
> WichHill(10)
 [1] 0.9168563 0.4666467 0.6701995 0.9274966 0.8431496 0.7882993 0.7696486
 [8] 0.2939858 0.6187752 0.9892436

in R:

> .Random.seed <- c(0, 1:3)  # select Wichmann-Hill, seed = 1:3
> runif(10)
 [1] 0.03381877 0.77754189 0.05273525 0.74462407 0.49036219 0.98285437
 [7] 0.80915099 0.71338138 0.80102091 0.98958603
> runif(10)
 [1] 0.9168563 0.4666467 0.6701995 0.9274966 0.8431496 0.7882993 0.7696486
 [8] 0.2939858 0.6187752 0.9892436


You can do this faster in C/Fortran, of course, but the load/save of the 
seed is slow from compiled code and in pure S I get over 1000 unifs/sec.

-- 
Brian D. Ripley,                  ripley at 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-help 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-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._