Générer des nombres aléatoires identiques en R et Julia
je voudrais générer des nombres aléatoires identiques en R et Julia. Les deux langues semblent utiliser la bibliothèque Mersenne-Twister par défaut, cependant dans Julia 1.0.0:
julia> using Random
julia> Random.seed!(3)
julia> rand()
0.8116984049958615
Produit 0.811..., tandis que dans R:
set.seed(3)
runif(1)
produit 0.168.
des idées?
Mon cas d'utilisation pour ceux qui sont intéressés: Tester de nouveaux Julia code nécessite la génération de nombres aléatoires (par exemple bootstrapping statistique) en comparant la sortie à celle des bibliothèques équivalentes dans R.
3 réponses
C'est un vieux problème.
Paul Gilbert a abordé la même question à la fin des années 1990 (!!) lorsqu'on essaie d'affirmer que les simulations en R (puis en newcomer) donnent le même résultat que celles en S-Plus (puis en Titulaire).
sa solution, et toujours L'approche dorée AFAICT: ré-implémenter en nouveau code dans les deux langues car c'est la seule façon d'assurer un semis identique, état, ... et tout ce qui l'affecte.
poursuivant le RCall suggestion faite par @Khashaa, il est clair que vous pouvez définir la graine et obtenir les nombres aléatoires de 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
et 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***
selon la suggestion de @ColinTBowers, voici une façon plus simple/plus propre d'y accéder R nombres aléatoires de 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
Voir:
?set.seed
"Mersenne-Twister": Tiré de Matsumoto et Nishimura (1998). Un GFSR torsadé avec période 2^19937 - 1 et répartition en 623 dimensions consécutives (sur l'ensemble de la période). Le 'seed' est un ensemble de 624 dimensions de 32 bits entiers plus une position courante dans cet ensemble.
et vous pourriez voir si vous pouvez faire un lien vers le même code C depuis les deux langues. Si vous voulez voir la liste/vecteur, tapez:
.Random.seed