来自随机正态分布的平行样本——不是更快吗?

Parallel samples from a random normal distribution - not faster?

我正在使用 R 创建一个从随机正态分布中获取样本的模拟,毫不奇怪,它相当慢。因此,我寻找一些使用 Rcpp 加速它的方法,并遇到了用于多线程计算的 RcppZiggurat package for faster random normal samples, and the RcppParallel package,我想,为什么不使用更快的算法并并行绘制样本呢?

于是开始做原型,最后总结出三种方法进行比较:

  1. 同时使用 RcppParallel 和 RcppZiggurat 的示例
  2. 仅使用 RcppZiggurat 的样本
  3. 样本使用旧 rnorm

下面是我使用 RcppParallel + RcppZiggurat(parallelDraws 函数)和 RcppZiggurat(serialDraws 函数)的实现:

#include <Rcpp.h>
// [[Rcpp::plugins("cpp11")]]
// [[Rcpp::depends(RcppParallel)]]
#include <RcppParallel.h>
// [[Rcpp::depends(RcppZiggurat)]]
#include <Ziggurat.h>

static Ziggurat::Ziggurat::Ziggurat zigg;

using namespace RcppParallel;

struct Norm : public Worker
{   
  int input;

  // saved draws
  RVector<double> draws;

  // constructors
  Norm(const int input, Rcpp::NumericVector draws)
    : input(input), draws(draws) {}

  void operator()(std::size_t begin, std::size_t end) {
    for (std::size_t i = begin; i < end; i++) {
      draws[i] = zigg.norm();
    }
  }
};

// [[Rcpp::export]]
Rcpp::NumericVector parallelDraws(int x) {

  // allocate the output vector
  Rcpp::NumericVector draws(x);

  // declare the Norm instance 
  Norm norm(x, draws);

  // call parallelFor to start the work
  parallelFor(0, x, norm);

  // return the draws
  return draws;
};

// [[Rcpp::export]]
Rcpp::NumericVector serialDraws(int x) {

  // allocate the output vector
  Rcpp::NumericVector draws(x);

  for (int i = 0; i < x; i++) {
    draws[i] = zigg.norm();
  }

  // return the draws
  return draws;
};

当我对它们进行基准测试时,我发现了一些令人惊讶的结果:

library(microbenchmark)
microbenchmark(parallelDraws(1e5), serialDraws(1e5), rnorm(1e5))
Unit: microseconds
                 expr      min       lq     mean    median       uq        max neval
 parallelDraws(1e+05) 3113.752 3539.686 3687.794 3599.1540 3943.282   5058.376   100
   serialDraws(1e+05)  695.501  734.593 2536.940  757.2325  806.135 175712.496   100
         rnorm(1e+05) 6072.043 6264.030 6655.835 6424.0195 6661.739  18578.669   100

单独使用 RcppZiggurat 比 rnorm 快大约 8 倍,但是一起使用 RcppParallel 和 RcppZiggurat 比单独使用 RcppZiggurat 慢!我尝试使用 RcppParallel ParallelFor 函数的 grain size,但它没有带来任何明显的改进。

我的问题是:在这里添加并行性实际上更糟的原因可能是什么?我知道 "overhead" 在并行计算中可能会超过各种因素带来的好处。这就是这里发生的事情吗?还是我完全误解了如何有效地使用 RcppParallel 包?

如评论中所述,开销可能会成问题,尤其是当整体运行时间较短时,最好不要将输出向量初始化为零并使用线程本地 RNG。示例实现:

#include <Rcpp.h>
// [[Rcpp::plugins("cpp11")]]
// [[Rcpp::depends(RcppParallel)]]
#include <RcppParallel.h>
// [[Rcpp::depends(RcppZiggurat)]]
#include <Ziggurat.h>


using namespace RcppParallel;

struct Norm : public Worker
{   
  // saved draws
  RVector<double> draws;

  // constructors
  Norm(Rcpp::NumericVector draws)
    : draws(draws) {}

  void operator()(std::size_t begin, std::size_t end) {
    Ziggurat::Ziggurat::Ziggurat zigg(end);
    for (std::size_t i = begin; i < end; i++) {
      draws[i] = zigg.norm();
    }
  }
};

// [[Rcpp::export]]
Rcpp::NumericVector parallelDraws(int x) {
  // allocate the output vector
  Rcpp::NumericVector draws(Rcpp::no_init(x));
  Norm norm(draws);
  parallelFor(0, x, norm);
  return draws;
}

// [[Rcpp::export]]
Rcpp::NumericVector serialDraws(int x) {
  // allocate the output vector
  Rcpp::NumericVector draws(Rcpp::no_init(x));
  Ziggurat::Ziggurat::Ziggurat zigg(42);
  for (int i = 0; i < x; i++) {
    draws[i] = zigg.norm();
  }
  return draws;
}

请注意,我使用的是 "poor man's parallel RNG",即不同线程的不同种子,希望一切顺利。我使用 end 作为种子,因为 begin 可能为零,我不确定 RcppZiggurat 中的 RNG 是否喜欢这样。由于创建 Ziggurat 对象需要一些时间(和内存),为了公平起见,我还使用本地对象进行串行计算。

对于 10^5 次随机抽取,使用并行计算仍然没有任何好处:

> bench::mark(parallelDraws(1e5), serialDraws(1e5), check = FALSE, min_iterations = 10)[,1:5]
# A tibble: 2 x 5
  expression                min   median `itr/sec` mem_alloc
  <bch:expr>           <bch:tm> <bch:tm>     <dbl> <bch:byt>
1 parallelDraws(1e+05)   1.08ms   1.78ms      558.     784KB
2 serialDraws(1e+05)   624.16µs  758.6µs     1315.     784KB

但是对于 10^8 次抽奖,我在我的双核笔记本电脑上得到了不错的 speed-up:

> bench::mark(parallelDraws(1e8), serialDraws(1e8), check = FALSE, min_iterations = 10)[,1:5]
# A tibble: 2 x 5
  expression                min   median `itr/sec` mem_alloc
  <bch:expr>           <bch:tm> <bch:tm>     <dbl> <bch:byt>
1 parallelDraws(1e+08)    326ms    343ms      2.91     763MB
2 serialDraws(1e+08)      757ms    770ms      1.30     763MB

因此,使用并行计算是否有意义在很大程度上取决于您需要的随机抽取次数。

顺便说一句,在我的 dqrng package is mentioned. This package also uses the Ziggurat method for normal (and exponential) draws combined with very fast 64bit RNGs, giving it comparable serial speed to RcppZiggurat for normal draws. In addition, the used RNGs are ment for parallel computation 评论中,即没有必要希望通过使用不同的种子来获得 non-overlapping 随机流。