Using rmultinom with Rcpp
bayesian, c++, mcmc, r, rcpp
Solution
Below are the answer by user95215 amended so that it compiles, and another version more in the Rcpp style:
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
IntegerVector oneMultinomC(NumericVector probs) {
int k = probs.size();
SEXP ans;
PROTECT(ans = Rf_allocVector(INTSXP, k));
probs = Rf_coerceVector(probs, REALSXP);
rmultinom(1, REAL(probs), k, &INTEGER(ans)[0]);
UNPROTECT(1);
return(ans);
}
// [[Rcpp::export]]
IntegerVector oneMultinomCalt(NumericVector probs) {
int k = probs.size();
IntegerVector ans(k);
rmultinom(1, probs.begin(), k, ans.begin());
return(ans);
}
Problem
I'd like to use the R function rmultinom in c++ code to be used with Rcpp. I get an error about not enough arguments - I am unfamiliar with what these arguments ought to be, as they do not corresond to the arguments used with the function in R. I also didn't have any luck using the "::Rf_foo" syntax to access an R function from Rcpp code. Here is a simplified version of my code below (yes, I am writing a gibbs sampler). ``` #include <Rcpp.h> using namespace Rcpp; // C++ implementation of the R which() function. int whichC(NumericVector x, double val) { int ind = -1; int n = x.size(); for (int i = 0; i < n; ++i) { if (x[i] == val) { if (ind == -1) { ind = i; } else { throw std::invalid_argument( "value appears multiple times." ); } } // end if } // end for if (ind != -1) { return ind; } else { throw std::invalid_argument( "value doesn't appear here!" ); return -1; } } // [[Rcpp::export]] int multSample(double p1, double p2, double p3) { NumericVector params(3); params(0) = p1; params(1) = p2; params(2) = p3; // HERE'S THE PROBLEM. RObject sampled = rmultinom(1, 1, params); int out = whichC(as<NumericVector>(sampled), 1); return out; } ``` I am new to c++, so I realize that alot of this code is probably nooby and inefficient. I'm open to suggestions on how to improve my c++ code, but my priority is hearing about the rmultinom business. Thanks! Btw, I apologize for similarity with this thread, but - The answer didn't work for my purposes - The difference may be enough to warrant a different question (do you think so? ) - The question was posted and answered a year ago.