[R] generating normal numbers: GetRNGstate, PutRNGstate

[Prof Brian Ripley]

On Wed, 3 Mar 2004, Yongchao Ge wrote:

> I'd like to generate thousands of normal numbers from my C function using 
> the C API functions provided R. I have two options:
> 
> 1. double norm_rand(); (page 61 of R extension 1.8.1)
> 2. double rnorm(double mu, double sigma); (page 58 of R extension 1.8.1)

Page numbers depend on the paper size you used (and I guess you didn't use
the world standard A4)  but what I guess is your page 61 says

  See Random numbers, for the protocol in using the random-variate 
  routines.

and that refers you back to page 58.

> If my understanding of R-exts is correct, then I only need to call 
> GetRNGstate once, and then call 1000 norm_rand, and then call 
> PutRNGstate once for the 1st option.
> 
> For the 2nd option, I have to call 1000 times for each of GetRNGstate, 
> rnorm, and PutRNGstate.

That is *not* what the manual says.  Just the same as 1).

> The pseudo-code for option 1 will be:
> 
> Method 1:
> 
> GetRNGstate();
> for(i=1;i<1000;i++){
>    x[i]=norm_rand();
> }
> PutRNGstate();
> 
> The pseudo code for option 2 will be:
> 
> Method 2:
> 
> for(i=1;i<1000;i++){
>    GetRNGstate();
>    x[i]=rnorm(0,1);
>    PutRNGstate();
> }
> 
> Of course, I can also write a slower version for option 1, i.e. call 
> GetRNGstate and PutRNGstate each time for norm_rand.
> 
> method 3:
> 
> for(i=1;i<1000;i++){
>    GetRNGstate();
>    x[i]=norm_rand(); 
>    PutRNGstate();
> }
> 
> 
> My questions are:
> 
> 1. Are the three methods all correct for generating random numbers?

Yes, given the answer to 2.

> 2. Are they generating the exactly the same random number if we have the 
> same random seed?

Yes.  Try it and see!

> I searched the R help and google and I didn't find answers. The reason for 
> to ask is that if both of the above answers are right, then I'd better 
> off use the method 1, which is the fastest as I have to generate 
> hundreds thousands of normal numbers.

See the above answer.  And please don't expect answers to questions like
this to be in the R-help: they are in the manuals and especially in the
source code.  It would be a good exercise for you to confirm this by

A) reading the R sources
B) doing some tests for yourself.

-- 
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 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595