CODEWITHSUNDEEP
rrcpp
Sacha Epskamp
Sacha Epskamp

Reputation: 47562

Use pnorm from Rmath.h with Rcpp

I am trying to write a piece of C++ code with Rcpp using functions like pnorm and qnorm. I can use the Rcpp sugar versions of these for vectors as explained in https://stackoverflow.com/a/9738848/567015, but I don't need to do this on a vector but just on a double.

If I understand correctly I can use the Rf_ prefix to get the scalar versions from Rmath.h. However, Rf_pnorm doesn't work:

library("inline")
Src <-  '
double y = as<double>(x);
double res = Rf_pnorm(y,0.0,1.0);
return wrap(res) ;
'

fx <- cxxfunction( signature(x = "numeric") ,Src, plugin = "Rcpp" )

fx(1)

Gives as errors:

file10c81a585dee.cpp: In function 'SEXPREC* file10c81a585dee(SEXP)':
file10c81a585dee.cpp:32:32: error: 'Rf_pnorm' was not declared in this scope

I found after some googling and trial and error that Rf_pnorm5 does work but requires extra parameters for lower tail and log scale:

Src <-  '
double y = as<double>(x);
double res = Rf_pnorm5(y,0.0,1.0,1,0);
return wrap(res) ;
'

fx <- cxxfunction( signature(x = "numeric") ,Src, plugin = "Rcpp" )

fx(1)
## [1] 0.8413447

Great, but I don't understand why this works but Rf_pnorm doesn't. I rather use Rf_pnorm because I think that makes it easier to find the right codes for different distributions.

Upvotes: 2

Views: 2290

Answers (1)

Dirk is no longer here
Dirk is no longer here

Reputation: 368261

Here is the Rcpp sugar variant which is more natural with Rcpp:

R> library(inline)
R> 
R> Src <- '
+ NumericVector y = NumericVector(x);
+ NumericVector res = pnorm(y,0.0,1.0);
+ return res;
+ '
R> 
R> fx <-  cxxfunction( signature(x = "numeric") , body=Src, plugin = "Rcpp")
R> 
R> fx(seq(0.8, 1.2, by=0.1))
[1] 0.788145 0.815940 0.841345 0.864334 0.884930
R> 
R> fx(1.0)      ## can also call with scalar args
[1] 0.841345
R> 

Looking more closely at our headers, we undefine the pnorm et al from Rmath.h in order to define the (vectorised) variants you get from Rcpp sugar.

Edit on 2012-11-14: With Rcpp 0.10.0 released today, you can call do the signature R::pnorm(double, double, double, int, int) if you want to use C-style code written against Rmath.h. Rcpp sugar still gives you vectorised versions.

Upvotes: 3

Related Questions