Home > database >  Poisson draw in Rcpp and R different results
Poisson draw in Rcpp and R different results

Time:01-29

I face the following contradiction when I use the same code in R and in Rcpp

In R I run the following code

t = 0
for(i in 1:50){
 t = t   rpois(1, 0.5)
}
t
[1] 28

and I take back a value t which is nonnegative. Now I type the exact same commands in Rcpp

#include <Rcpp.h>
#include<Rmath.h>
using namespace Rcpp;

// [[Rcpp::export]]
int Pois(int l){ 
  int t=0;
  for(int i=0; i<50;  i){
    t =R::rpois(l);
  }
  return t;
}

and when I call the function in R

Pois(0.5)
[1] 0

which is wrong since in R it was different of zero

What is going wrong?

CodePudding user response:

You should use double l rather than int l, e.g.,

int Pois(double l){ 
  int t=0;
  for(int i=0; i<50;  i){
    t =R::rpois(l);
  }
  return t;
}

otherwise (int) 0.5 gives you 0.

CodePudding user response:

@ThomasIsCoding already showed you the main issue. But recall that beside R::rpois() we also have the vectorised Rcpp::rpois(). And, as usual, given the same seed it gives the same draws as R:

> set.seed(123)
> rpois(10, 0.5)
 [1] 0 1 0 1 2 0 0 1 0 0
> Rcpp::cppFunction("NumericVector myrp(int n, double l) { return Rcpp::rpois(n, l); }")
> set.seed(123)
> myrp(10, 0.5)
 [1] 0 1 0 1 2 0 0 1 0 0
> 
  •  Tags:  
  • Related