The function runif(n) contains a protection against a corrupted .Random.seed.
Besides other things, it tests whether at least one of the numbers
.Random.seed[3:626] is not zero. If all(.Random.seed[3:626]==0), then
the internal Mersenne Twister state is regenerated using current time,
since a zero state is a fixed point of the recurrence and, hence, produces
a constant sequence.
However, the condition any(.Random.seed[3:626]!=0) does not imply that
the internal Mersenne Twister state is indeed not zero, since only the
most significant bit of .Random.seed[3] belongs to the internal state.
Hence, the number of bits in the state of Mersenne Twister
is 624*32 - 31 = 19937, which explains the period 2^19937-1.
For example, if .Random.seed[3] == 1, we always have the condition
any(.Random.seed[3:626]!=0)
satisfied, but the internal state may still be effectively zero. An
example of such a situation is
RNGkind("default")
.Random.seed[3:626] <- as.integer(0)
.Random.seed[3] <- as.integer(1)
x <- runif(10000)
all(x==x[1]) # TRUE
length(unique(x)) # 1
all(.Random.seed[3:626]==0) # TRUE
Here, the internal state was effectively zero, but this fact was not
detected, since some (unimportant) bits of .Random.seed[3] were not zero.
On the contrary, if also .Random.seed[3]==0, then the internal state
is regenerated and the output of runif() becomes non constant:
RNGkind("default")
.Random.seed[3:626] <- as.integer(0)
x <- runif(10000)
all(x==x[1]) # FALSE
length(unique(x)) # 10000
all(.Random.seed[3:626]==0) # FALSE
The following patch to FixupSeeds corrects the detection of zero state:
--- R-devel_2007-10-14-orig/src/main/RNG.c 2007-09-02 07:49:35.000000000
+0200
+++ R-devel_2007-10-14-FixupSeeds/src/main/RNG.c 2007-10-15
04:33:52.988060624 +0200
@@ -181,7 +181,9 @@
/* No action unless user has corrupted .Random.seed */
if(I1 <= 0) I1 = 624;
/* check for all zeroes */
- for (j = 1; j <= 624; j++)
+ notallzero = ((RNG_Table[RNG_kind].i_seed[1] & 0x80000000) != 0);
+ if(!notallzero)
+ for (j = 2; j <= 624; j++)
if(RNG_Table[RNG_kind].i_seed[j] != 0) {
notallzero = 1;
break;
Petr Savicky.
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
