I'd like to generate identical random numbers in R and Julia. Both languages appear to use the Mersenne-Twister library by default, however in Julia 1.0.0:
julia> using Random
julia> Random.seed!(3)
julia> rand()
0.8116984049958615
Produces 0.811...
, while in R:
set.seed(3)
runif(1)
produces 0.168
.
Any ideas?
Related SO questions here and here.
My use case for those who are interested: Testing new Julia code that requires random number generation (e.g. statistical bootstrapping) by comparing output to that from equivalent libraries in R.
In Julia, you can set a seed to the random number generator using the srand() function. The code example below sets the seed to 1234. Generating a random variable with rand(1) after setting the seed to 1234 will always generate the same number, i.e. it will always return 0.5908446386657102.
Use rand(Uniform(-2,1),3,3) or rand(Normal(x,y), 3,3) to get normally distributed values.
That is an old problem.
Paul Gilbert addressed the same issue in the late 1990s (!!) when trying to assert that simulations in R (then then newcomer) gave the same result as those in S-Plus (then the incumbent).
His solution, and still the golden approach AFAICT: re-implement in fresh code in both languages as the this the only way to ensure identical seeding, state, ... and whatever else affects it.
Pursuing the RCall
suggestion made by @Khashaa, it's clear that you can set the seed and get the random numbers from R
.
julia> using RCall
julia> RCall.reval("set.seed(3)")
RCall.NilSxp(16777344,Ptr{Void} @0x0a4b6330)
julia> a = zeros(Float64,20);
julia> unsafe_copy!(pointer(a), RCall.reval("runif(20)").pv, 20)
Ptr{Float64} @0x972f4860
julia> map(x -> @printf("%20.15f\n", x), a);
0.168041526339948
0.807516399072483
0.384942351374775
0.327734317164868
0.602100674761459
0.604394054040313
0.124633444240317
0.294600924244151
0.577609919011593
0.630979274399579
0.512015897547826
0.505023914156482
0.534035353455693
0.557249435689300
0.867919487645850
0.829708693316206
0.111449153395370
0.703688358888030
0.897488264366984
0.279732553754002
and from R
:
> options(digits=15)
> set.seed(3)
> runif(20)
[1] 0.168041526339948 0.807516399072483 0.384942351374775 0.327734317164868
[5] 0.602100674761459 0.604394054040313 0.124633444240317 0.294600924244151
[9] 0.577609919011593 0.630979274399579 0.512015897547826 0.505023914156482
[13] 0.534035353455693 0.557249435689300 0.867919487645850 0.829708693316206
[17] 0.111449153395370 0.703688358888030 0.897488264366984 0.279732553754002
** EDIT **
Per the suggestion by @ColinTBowers, here's a simpler/cleaner way to access R
random numbers from Julia
.
julia> using RCall
julia> reval("set.seed(3)");
julia> a = rcopy("runif(20)");
julia> map(x -> @printf("%20.15f\n", x), a);
0.168041526339948
0.807516399072483
0.384942351374775
0.327734317164868
0.602100674761459
0.604394054040313
0.124633444240317
0.294600924244151
0.577609919011593
0.630979274399579
0.512015897547826
0.505023914156482
0.534035353455693
0.557249435689300
0.867919487645850
0.829708693316206
0.111449153395370
0.703688358888030
0.897488264366984
0.279732553754002
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With