Функция R replicate () не работает с функцией Rcpp - PullRequest
Новая
       11

Функция R replicate () не работает с функцией Rcpp

0 голосов
/ 26 мая 2018

У меня проблемы с использованием функции replicate() в R для генерации случайных чисел с функцией Rcpp.Рассмотрим следующую функцию в R: 

trial <- function(){rnorm(1)}
replicate(10, trial())

. Она генерирует 10 случайных чисел из распределения Гаусса.Он работает совершенно нормально и выдает результат, подобный следующему: 

 [1]  0.7609912 -0.2949613  1.8684363 -0.3358377 -1.6043926  0.2706250  0.5528813  1.0228125 -0.2419092 -1.4761937

Однако у меня есть функция c ++ getRan(), которая генерирует случайное число из распределения Гаусса.Я снова использовал replicate для вызова функции следующим образом: 

replicate(10,getRan())

Он создает вектор с таким же номером, как этот: 

> replicate(10,getRan())
 [1] -1.374932 -1.374932 -1.374932 -1.374932 -1.374932 -1.374932 -1.374932 -1.374932 -1.374932 -1.374932
> replicate(10,getRan())
 [1] -0.3273785 -0.3273785 -0.3273785 -0.3273785 -0.3273785 -0.3273785 -0.3273785 -0.3273785 -0.3273785 -0.3273785
> replicate(10,getRan())
 [1] -0.7591953 -0.7591953 -0.7591953 -0.7591953 -0.7591953 -0.7591953 -0.7591953 -0.7591953 -0.7591953 -0.7591953
> replicate(10,getRan())
 [1] -1.698935 -1.698935 -1.698935 -1.698935 -1.698935 -1.698935 -1.698935 -1.698935 -1.698935 -1.698935

Однако, если я вызываю функцию несколько раз, онаотлично работает: 

 getRan()
[1] 1.345227
> getRan()
[1] 0.3555393
> getRan()
[1] 1.587241
> getRan()
[1] 0.5313518

Так в чем здесь проблема?Повторяет ли функция replicate() одну и ту же функцию, возвращаемую из getRan(), вместо вызова getRan() несколько раз?Это ошибка?

PS: я знаю, что могу использовать rnorm(n) для генерации n нормальных случайных чисел, но я хочу использовать функцию c ++ для более сложных вычислений, основанных на генерации случайных чисел

PPS: это мой код на C ++: 

double getRan(){
  unsigned seed = std::chrono::system_clock::now().time_since_epoch().count();
  std::default_random_engine generator(seed);
  std::normal_distribution<double> distribution (0.0,1.0);
  double epi = distribution(generator);
  return epi;
}

Ответы [ 2 ]

0 голосов
/ 29 мая 2018

Диркс ответ правильный.Вы должны использовать R's RNG.Если вы настаиваете на использовании RNG в C ++, вы можете использовать что-то вроде этого: 

#include <Rcpp.h>
// [[Rcpp::plugins(cpp11)]]
#include <random>

namespace {
  std::default_random_engine generator(std::random_device{}());
  std::normal_distribution<double> distribution (0.0,1.0);  
}

// [[Rcpp::export]]
double getRan(){
  return distribution(generator);
}

/*** R
replicate(10,getRan())
*/

Это позволяет избежать создания нового экземпляра std::default_random_engine (и std::normal_distribution) при каждом вызове функции.Это важно, поскольку свойства ГСЧ гарантируются только для повторных отборов из одного ГСЧ.Не для повторных розыгрышей из разных ГСЧ, засеянных (возможно, разными) семенами.

Кстати, в моей системе ваш исходный код не выдает одно и то же число несколько раз.Если у вас проблемы с std::random_device и вы работаете с Windows, вы можете быть подвержены этой ошибке mingw .В этом случае семена по времени - лучшая альтернатива.

0 голосов
/ 26 мая 2018

Вот контрпример, показывающий, что он работает просто отлично:

Код 

trialR <- function() { rnorm(1) }
Rcpp::cppFunction("double trialC() { return R::rnorm(0.0, 1.0); }")
Rcpp::cppFunction("Rcpp::NumericVector trialSugar() { return Rcpp::rnorm(1.0, 0.0, 1.0); }")

set.seed(123); replicate(3, trialR())
set.seed(123); replicate(3, trialC())
set.seed(123); replicate(3, trialSugar())

Выход

Через Rscript для обеспечения нового сеанса и т. Д. Pp:

edd@rob:/tmp$ Rscript so50543659.R 
[1] -0.560476 -0.230177  1.558708
[1] -0.560476 -0.230177  1.558708
[1] -0.560476 -0.230177  1.558708
edd@rob:/tmp$
...