On 29 October 2014 at 14:48, Helske Jouni wrote: | Dear list, | | | | I am trying to call R’s negative binomial function with Rcpp, but encountered | some weird behaviour. Here is an example using inline: | | | | | rcpp_Rf_dnbinom <- rcpp(signature(), | ' return wrap(Rf_dnbinom( 4.0, 0.5, 0.9, 1)); ') | | rcpp_Rf_dnbinom_mu <- rcpp(signature(), | ' return wrap(Rf_dnbinom_mu( 4.0, 0.5, 0.9, 1)); ') | | | rcpp_dnbinom <- rcpp(signature(), | ' return wrap(R::dnbinom( 4.0, 0.5, 0.9, 1)); ') | | rcpp_dnbinom_mu <- rcpp(signature(), | ' return wrap(R::dnbinom_mu( 4.0, 0.5, 0.9, 1)); ') | | rcpp_dnbinom_sugar <- rcpp(signature(y="numeric"), | ' NumericVector x = NumericVector(y); | NumericVector res = dnbinom( x, 0.5, 0.9,1); | return wrap(res); ') | | | rcpp_dnbinom_mu_sugar <- rcpp(signature(y="numeric"), | ' NumericVector x = NumericVector(y); | NumericVector res = dnbinom_mu( x, 0.5, 0.9,1); | return wrap(res); ') | | rcpp_Rf_dnbinom() | #-10.5597 | rcpp_Rf_dnbinom_mu() | #-3.578823 | rcpp_dnbinom() | #-10.5597 | rcpp_dnbinom_mu() | #-10.5597 | rcpp_dnbinom_sugar(y=4.0) | #-10.5597 | rcpp_dnbinom_mu_sugar(y=4.0) | #-3.578823 | dnbinom(x=4,size=0.5,mu=0.9,log=TRUE) | #[1] -3.578823 | dnbinom(x=4,size=0.5,prob=0.9,log=TRUE) | #[1] -10.5597 | | | | So it looks like that everything is fine when using Rcpp sugar or Rf_dbinom_mu | directly, but when using form R::dbinom both dnbinom and dnbinom_mu calls | actually use dnbinom. | | | I am not sure if this is relevant, but Rmath.h in Rcpp contains lines: | | | | 129 inline double dnbinom(double x, double sz, double pb, int lg) { return :: | Rf_dnbinom(x, sz, pb, lg); } | | 130 inline double pnbinom(double x, double sz, double pb, int lt, int lg) { | return ::Rf_pnbinom(x, sz, pb, lt, lg); } | | 131 inline double qnbinom(double p, double sz, double pb, int lt, int lg) { | return ::Rf_qnbinom(p, sz, pb, lt, lg); } | | 132 inline double rnbinom(double sz, double pb) { return ::Rf_rnbinom(sz, pb); | } | | 133 | | 134 inline double dnbinom_mu(double x, double sz, double mu, int lg) { return | ::Rf_dnbinom(x, sz, mu, lg); } | | 135 inline double pnbinom_mu(double x, double sz, double mu, int lt, int lg) { | return ::Rf_pnbinom(x, sz, mu, lt, lg); } | | 136 inline double qnbinom_mu(double x, double sz, double mu, int lt, int lg) { | return ::Rf_qnbinom(x, sz, mu, lt, lg); } | | 137 inline double rnbinom_mu(double sz, double mu) { return ::Rf_rnbinom(sz, | mu); } | | | | To me this looks like that both dnbinom and dnbinom_mu calls the same | function.
That is likely a bug due an oversight of mine. The four lines in 134 to 137 want to call the variants from R ending in _mu. | So, is there something wrong with my code or is there a typo in R:: dnbinom_mu? | | | To be honest, I don't really understand the differences between the three See eg help(qbinom): one can either supply mu or prob when using the _R_ variant, so when I wrote this interface I tried to mimic this, offer the differnt functions with and without mu (in which case prob is used) but then dropped the ball amd didn't add _mu in the call. Will fix -- t hanke for the bug report! | versions, except that the sugar version needs NumericVector as first argument | and it is vectorized regards that parameter. I would actually need version | which is vectorized wrt all arguments (except log), but since there isn't one I | am doing something like this: | | | #include "RcppArmadillo.h" | | // [[Rcpp::depends(RcppArmadillo)]] | | double fun(const arma::mat& y, const arma::mat& x ,const arma::mat& theta){ | double res=0.0; | for(unsigned int i=0; i<y.n_elem; i++){ | if(arma::is_finite(y(i))){ | res += Rf_dbinom_mu( y(i), u(i), theta(i), 1); //was R::dbinom_mu | } | } | | return res; | | } I write similar little helpers when needed. Our entire sugar interface generally has vector 'x' but scalar auxiliary parameters. Dirk | | Best regards, | | | | Jouni Helske | | | | _______________________________________________ | Rcpp-devel mailing list | Rcpp-devel@lists.r-forge.r-project.org | https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel -- http://dirk.eddelbuettel.com | @eddelbuettel | e...@debian.org _______________________________________________ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel