Skip to content

[Rcpp-devel] R random numbers vs c++

6 messages · Glenn Lawyer, Dirk Eddelbuettel, Davor Cubranic +1 more

#
Is there a reason to prefer Rcpp::runif over the c++ rand() for, 
say, accepting a state in a Markov chain simulation? Is one more 
random than another? faster?

+glenn
#
On 22 March 2012 at 13:41, Glenn Lawyer wrote:
| Is there a reason to prefer Rcpp::runif over the c++ rand() for, say, accepting
| a state in a Markov chain simulation? Is one more random than another? faster?

Fair question.  I'd say go with R because

 - it is portable as R supplies it, and well tested by R Core (specifically
   Brian Ripley who actually "did write the book on this")

 - six different generators you can select from

 - numerous distributions

whereas system library rand() implementations used to be terrible.  But you
could also consider Boost or others.

As for speed, we were benchmarking some MCMC stuff against Python in response
to a post by Darren Wilkinson last year. Turned out that R's rgamma was
slower that one from the GSL, which I discussed briefly with Martin Maechler.
It "may" be that R's rgamma is more rigorous / better.

Dirk
#
On 2012-03-22, at 5:41 AM, Glenn Lawyer wrote:

            
I'd use the R random generator so your C++ code is in sync with R in terms of seeds and random sequences. For instance, if you're running your simulation in parallel on multiple machines or cores, you can use one of the RNG's designed specifically to avoid correlations between individual instances, such as SPRNG or RNGstream.

Davor
#
On 22 March 2012 at 08:20, Davor Cubranic wrote:
| On 2012-03-22, at 5:41 AM, Glenn Lawyer wrote:
| 
| > Is there a reason to prefer Rcpp::runif over the c++ rand() for, say, accepting a state in a Markov chain simulation? Is one more random than another? faster?
| 
| I'd use the R random generator so your C++ code is in sync with R in terms of seeds and random sequences. For instance, if you're running your simulation in parallel on multiple machines or cores, you can use one of the RNG's designed specifically to avoid correlations between individual instances, such as SPRNG or RNGstream.

Ah, yes, indeed. Sorry I omitted that. It is indeed important that

    set.seed(42); rnorm(6)

and

    f <- cxxfunction(,body='RNGScope tmp; return rnorm(2)', plugin="Rcpp")
    set.seed(42); c(rnorm(2), f(), rnorm(2))

return the same result:

    R> f <- inline::cxxfunction(,body='RNGScope tmp; return rnorm(2);', plugin="Rcpp")
    R> set.seed(42); c(rnorm(2), f(), rnorm(2))
    [1]  1.370958 -0.564698  0.363128  0.632863  0.404268 -0.106125
    R> 
    R> set.seed(42); rnorm(6)
    [1]  1.370958 -0.564698  0.363128  0.632863  0.404268 -0.106125
    R> 

so that you can mix and match.

Dirk
#
FWIW, the random-number stuff in C++11 appears to be greatly expanded:

??? http://en.wikipedia.org/wiki/C%2B%2B11#Extensible_random_number_facility

This is discussed in some more detail in, for instance:

??? http://www.amazon.com/Professional-C-Wrox-Guides/dp/0470932449/ref=sr_1_1?ie=UTF8&qid=1332436173&sr=8-1

It's my understanding the Rcpp does not (yet?) support C++11.

-- Mike
#
On 22 March 2012 at 10:13, Michael Hannon wrote:
| >Is there a reason to prefer Rcpp::runif over the c++ rand() for, say,
| >accepting a state in a Markov chain simulation? Is one more random than
| >another? faster?
| 
| FWIW, the random-number stuff in C++11 appears to be greatly expanded:
| 
| ??? http://en.wikipedia.org/wiki/C%2B%2B11#Extensible_random_number_facility
| 
| This is discussed in some more detail in, for instance:
| 
| ??? http://www.amazon.com/Professional-C-Wrox-Guides/dp/0470932449/ref=sr_1_1?ie=UTF8&qid=1332436173&sr=8-1
| 
| It's my understanding the Rcpp does not (yet?) support C++11.

It is a bit more nuanced. 

Romain and I jumped on a few C++11 features early but then realized that we
are bound by what CRAN has as its standard compiler setup.  And that does not
yet support C++11.  But you can of course use it locally if you add the g++
flag "-std=c++0x" --- simply don't upload a src/Makevars with that to CRAN.

With this switch, you then get this activated in RcppCommon.h:

#ifdef __GNUC__
    #define GCC_VERSION (__GNUC__ * 10000 + __GNUC_MINOR__ * 100 + __GNUC_PATCHLEVEL__)
    #ifdef __GXX_EXPERIMENTAL_CXX0X__
        #define HAS_CXX0X
        #if GCC_VERSION >= 40300
            #define HAS_VARIADIC_TEMPLATES
        #endif
        #if GCC_VERSION >= 40400
            #define HAS_INIT_LISTS
        #endif
    #endif

so if a recent g++ with -std=c++0x (which gets us __GXX_EXPERIMENTAL_CXXOX_)
is used, HAS_CXX0X is active.  

You can grep around and see what that gets tested for that.  We haven't added
much recently due to the "hold" on C++0x at CRAN.  

But to come back to the thread, yes lots of nice new things in C++11. Not
your grandma's C++.  Now the compilers just need to catch up...

Dirk