[R] difference between rnorm(1000, 0, 1) and running rnorm(500, 0, 1) twice

[Duncan Murdoch]

On 2/8/2006 8:30 AM, Brian D Ripley wrote:
> On Wed, 8 Feb 2006, Duncan Murdoch wrote:
> 
>> On 2/8/2006 4:53 AM, Bjï¿½rn-Helge Mevik wrote:
>> > Why don't you test it yourself?
>> >
>> > E.g.,
>> >
>> > set.seed(42)
>> > bob1 <- rnorm(1000,0,1)
>> > set.seed(42)
>> > bob2 <- rnorm(500,0,1)
>> > bob3 <- rnorm(500,0,1)
>> > identical(bob1, c(bob2, bob3))
>> >
>> > I won't tell you the answer. :-)
>>
>> This isn't really something that can be proved by a test.  Perhaps the
>> current implementation makes those equal only because 500 is even, or
>> divisible by 5, or whatever...
>>
>> I think the intention is that those should be equal, but in a quick
>> search I've been unable to find a documented guarantee of that.  So I
>> would take a defensive stance and assume that there may be conditions
>> where c(rnorm(m), rnorm(n)) is not equal to rnorm(m+n).
>>
>> If someone can point out the document I missed, I'd appreciate it.
> 
> It's various source files in R_HOME/src/main.
> 
> Barring bugs, they will be the same.  As you know
> 
> 	R is free software and comes with ABSOLUTELY NO WARRANTY.

I didn't mean guarantee in the sense of warranty, just guarantee in the 
sense that if someone found a situation where they weren't equal, we 
would consider it a bug and fix it or document it as an exception.

Should we add a statement to the RNG man page or manuals somewhere that 
says this is the intention?

For others who aren't as familiar with the issues as Brian: this isn't 
necessarily a good idea.  We have a lot of RNGs, and it's fairly easy to 
write one so that this isn't true.  For example, the Box-Muller method 
naturally generates pairs of normals; a naive implementation would just 
throw one away at the end if asked for an odd number.  (Ours doesn't do 
that.)

Duncan Murdoch