RcppArmadillo 伽玛分布在具有相同种子的平台之间不同
RcppArmadillo gamma distribution differs between platforms with same seed
我正在研究 a package,它使用来自 RcppArmadillo 的随机数。该软件包运行 MCMC 算法,为了获得精确的再现性,用户应该能够设置随机数种子。执行此操作时,似乎 arma::randg()
函数用于从伽马分布 returns 跨平台的不同值生成随机数。 arma::randu()
或 arma::randn()
不是这种情况。这可能与 arma::randg()
需要 C++11 这一事实有关吗?
这是我在 macOS 10.13.6、运行 R3.5.2 上得到的:
library(Rcpp)
library(RcppArmadillo)
sourceCpp(code = '
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
double random_gamma() {
return arma::randg();
}
// [[Rcpp::export]]
double random_uniform() {
return arma::randu();
}
// [[Rcpp::export]]
double random_normal() {
return arma::randn();
}
'
)
replicate(2, {set.seed(1); random_gamma()})
#> [1] 1.507675 1.507675
replicate(2, {set.seed(432); random_gamma()})
#> [1] 0.02234341 0.02234341
replicate(2, {set.seed(1); random_uniform()})
#> [1] 0.2655087 0.2655087
replicate(2, {set.seed(1); random_normal()})
#> [1] -1.390378 -1.390378
由 reprex package (v0.2.1)
于 2019-02-22 创建
这是我在 Windows 10, 运行 R3.5.2:
上得到的
library(Rcpp)
library(RcppArmadillo)
sourceCpp(code = '
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
double random_gamma() {
return arma::randg();
}
// [[Rcpp::export]]
double random_uniform() {
return arma::randu();
}
// [[Rcpp::export]]
double random_normal() {
return arma::randn();
}
'
)
replicate(2, {set.seed(1); random_gamma()})
#> [1] 0.2549381 0.2549381
replicate(2, {set.seed(432); random_gamma()})
#> [1] 0.2648896 0.2648896
replicate(2, {set.seed(1); random_uniform()})
#> [1] 0.2655087 0.2655087
replicate(2, {set.seed(1); random_normal()})
#> [1] -1.390378 -1.390378
由 reprex package (v0.2.1)
于 2019-02-22 创建
可以看出,arma::randg()
生成的随机数内部是一致的,但平台不同。
我尝试使用 Armadillo 文档中的说明设置种子:
library(Rcpp)
library(RcppArmadillo)
sourceCpp(code = '
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
double random_gamma(int seed) {
arma::arma_rng::set_seed(seed);
return arma::randg();
}
'
)
replicate(4, random_gamma(1))
#> Warning in random_gamma(1): When called from R, the RNG seed has to be set
#> at the R level via set.seed()
#> [1] 1.3659195 0.6447221 1.1771862 0.9099034
由 reprex package (v0.2.1)
于 2019-02-22 创建
但是,如警告和结果所示,种子不会以这种方式设置。
有没有办法在使用 arma::randg()
时在平台之间获得可重现的结果,或者我是否需要使用 RcppArmadillo 中可用的其他一些随机数生成器来实现伽玛分布?
更新
正如评论中指出的那样,使用 R::rgamma()
解决了这个问题。以下代码 returns Mac 和 Windows 上的相同数字:
library(Rcpp)
sourceCpp(code = '
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
double random_gamma() {
return R::rgamma(1.0, 1.0);
}
'
)
replicate(2, {set.seed(1); random_gamma()})
#> [1] 0.1551414 0.1551414
由 reprex package (v0.2.1)
于 2019-02-22 创建
这帮我解决了。但是,我不确定问题是否已解决,因为这似乎是意外行为,所以保持打开状态。
总结评论中的讨论:
- 对于 gamma 分布,Armadillo 使用 C++11 的
random
header、c.f 中的 std::gamma_distribution
。 https://gitlab.com/conradsnicta/armadillo-code/blob/9.300.x/include/armadillo_bits/arma_rng_cxx11.hpp#L165
- 在 C++ 中生成标准随机数分布的算法是 implementation-defined。
- 如果您需要 cross-platform 可重现的结果,最简单的解决方案是使用通过
R::rgamma
或 Rcpp::rgamma
在 R 中实现的伽玛分布。
我正在研究 a package,它使用来自 RcppArmadillo 的随机数。该软件包运行 MCMC 算法,为了获得精确的再现性,用户应该能够设置随机数种子。执行此操作时,似乎 arma::randg()
函数用于从伽马分布 returns 跨平台的不同值生成随机数。 arma::randu()
或 arma::randn()
不是这种情况。这可能与 arma::randg()
需要 C++11 这一事实有关吗?
这是我在 macOS 10.13.6、运行 R3.5.2 上得到的:
library(Rcpp)
library(RcppArmadillo)
sourceCpp(code = '
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
double random_gamma() {
return arma::randg();
}
// [[Rcpp::export]]
double random_uniform() {
return arma::randu();
}
// [[Rcpp::export]]
double random_normal() {
return arma::randn();
}
'
)
replicate(2, {set.seed(1); random_gamma()})
#> [1] 1.507675 1.507675
replicate(2, {set.seed(432); random_gamma()})
#> [1] 0.02234341 0.02234341
replicate(2, {set.seed(1); random_uniform()})
#> [1] 0.2655087 0.2655087
replicate(2, {set.seed(1); random_normal()})
#> [1] -1.390378 -1.390378
由 reprex package (v0.2.1)
于 2019-02-22 创建这是我在 Windows 10, 运行 R3.5.2:
上得到的library(Rcpp)
library(RcppArmadillo)
sourceCpp(code = '
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
double random_gamma() {
return arma::randg();
}
// [[Rcpp::export]]
double random_uniform() {
return arma::randu();
}
// [[Rcpp::export]]
double random_normal() {
return arma::randn();
}
'
)
replicate(2, {set.seed(1); random_gamma()})
#> [1] 0.2549381 0.2549381
replicate(2, {set.seed(432); random_gamma()})
#> [1] 0.2648896 0.2648896
replicate(2, {set.seed(1); random_uniform()})
#> [1] 0.2655087 0.2655087
replicate(2, {set.seed(1); random_normal()})
#> [1] -1.390378 -1.390378
由 reprex package (v0.2.1)
于 2019-02-22 创建可以看出,arma::randg()
生成的随机数内部是一致的,但平台不同。
我尝试使用 Armadillo 文档中的说明设置种子:
library(Rcpp)
library(RcppArmadillo)
sourceCpp(code = '
#include <RcppArmadillo.h>
using namespace Rcpp;
// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
double random_gamma(int seed) {
arma::arma_rng::set_seed(seed);
return arma::randg();
}
'
)
replicate(4, random_gamma(1))
#> Warning in random_gamma(1): When called from R, the RNG seed has to be set
#> at the R level via set.seed()
#> [1] 1.3659195 0.6447221 1.1771862 0.9099034
由 reprex package (v0.2.1)
于 2019-02-22 创建但是,如警告和结果所示,种子不会以这种方式设置。
有没有办法在使用 arma::randg()
时在平台之间获得可重现的结果,或者我是否需要使用 RcppArmadillo 中可用的其他一些随机数生成器来实现伽玛分布?
更新
正如评论中指出的那样,使用 R::rgamma()
解决了这个问题。以下代码 returns Mac 和 Windows 上的相同数字:
library(Rcpp)
sourceCpp(code = '
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
double random_gamma() {
return R::rgamma(1.0, 1.0);
}
'
)
replicate(2, {set.seed(1); random_gamma()})
#> [1] 0.1551414 0.1551414
由 reprex package (v0.2.1)
于 2019-02-22 创建这帮我解决了。但是,我不确定问题是否已解决,因为这似乎是意外行为,所以保持打开状态。
总结评论中的讨论:
- 对于 gamma 分布,Armadillo 使用 C++11 的
random
header、c.f 中的std::gamma_distribution
。 https://gitlab.com/conradsnicta/armadillo-code/blob/9.300.x/include/armadillo_bits/arma_rng_cxx11.hpp#L165 - 在 C++ 中生成标准随机数分布的算法是 implementation-defined。
- 如果您需要 cross-platform 可重现的结果,最简单的解决方案是使用通过
R::rgamma
或Rcpp::rgamma
在 R 中实现的伽玛分布。