Transfert de semences/Définition de semences/C dans le code R

Question

Je fais C dans les codes R.

Dans mon code C, j'utilise la fonction rand() pour générer un nombre aléatoire. Le R-ext.pdf dit que je dois définir une graine en utilisant les commandes; Bien que j'utilise ces commandes ci-dessus, je reçois toujours des valeurs différentes pour la même graine. Pourriez-vous me donner de l'aide?

L'exemple minimum est:

En C:

# include <R.h> 
# include <Rinternals.h> 
# include <Rmath.h> 
# include <R_ext/Linpack.h> 

SEXP example(){ 

    SEXP output; 
    GetRNGstate(); 
    PROTECT(output = allocVector(INTSXP, 1)); 
    INTEGER(output)[0] = rand() % 50; 
    PutRNGstate(); 
    UNPROTECT(1); 
    return(output); 
} 

En R:

dyn.load("example.so") 
## The following codes return different values at ever run 
set.seed(1) 
.Call("example") 

Merci à l'avance.

Answer 1

C'est une erreur de logique dans votre pensée - vous correctement définir la graine, initialiser le R RNG à partir du code ... mais appelez le RNG système au lieu de la R RNG.

Remplacez rand() par unif_rand() (ou norm_rand()) et vous devriez être activé.

Rcpp rend tout cela plus facile et vous donne un accès vectorisé aux tirages de la fonction de distribution (mais vous pouvez bien sûr faire tout cela à la main en C si vous préférez).

En utilisant cppFunction() de Rcpp, nous maintenant prendre soin de RNGScope qui fournit à son tour GetRNGstate()/PutRNGstate() (alors que des exemples plus anciens montrent encore instanciation de RNGScope, ajoutant ne fait aucun mal comme il le fait un équivalent de comptage de référence).

Il est vraiment définir une doublure, étendre automatiquement, compiler et charger ceci:

R> cppFunction("double myrand() { return norm_rand(); }") 
R> for (i in 1:5) { set.seed(42); cat(i, " -- ", myrand(), "\n") } 
1 -- 1.37096 
2 -- 1.37096 
3 -- 1.37096 
4 -- 1.37096 
5 -- 1.37096 
R> 

alors sans réensemencement nous obtenons

R> for (i in 1:5) { cat(i, " -- ", myrand(), "\n") } 
1 -- -0.564698 
2 -- 0.363128 
3 -- 0.632863 
4 -- 0.404268 
5 -- -0.106125 
R> 

Enfin, si vous vraiment voulez vous pouvez bien sûr continuer à utiliser rand() (mais voir la littérature sur ses performances sucky), mais utilisez sa fonction d'ensemencement au lieu de R.