[Rcpp-devel] R.e. Fwd: problem with rmultinom function
On Tue, May 31, 2011 at 10:55 AM,
<rcpp-devel-request at r-forge.wu-wien.ac.at> wrote:
Hi, I sent this message a couple of days ago, but I haven't heard anything yet. I think you didnt receive it, that's why I send it again.
I didn't see anything. If you're not subscribed to the list, submissions silently fail. Plus, you must send from the subscribing address. I've made this mistake a few times...
I have a strange problem that I spent several days trying to fix it, but unfortunately I wasn't able to catch it. I'm running the rcpp code for 3000 iterations(running the same code for 3000 times) on different datasets. but for some datasets, in the final iterations (like 2667 th iteration) the code stops running and it jumps out of the R environment. and if I run the same code (without any changes) on the same dataset(still without any changes) the problem might never happen and it goes to the last iteration without any problem(or the problem may occur in another iteration). and the interesting part is that error code is different each time, I got error codes like: what(): ?invalid partial string match what(): ?'prob' is missing what(): ?unused argument(s) (<environment>) and etc... I'm pretty sure that problem occurs in the line of code that I am using rmultinom function, but I included the necessary function Function rmultinom = stats["rmultinom"]; to be more clear I wrote a simple function and I attached it in this email. in this code I didn't read any dataset, I just used some variables to define the size of structures that I need to get from dataset. if you run the code with these values for P, K, unum,qnum, the error occurs in the first 10 iterations, but If you change the values of those parameters to smaller values, the error may never happen. I tried to used an cpp implementation for rmultinom, it works but the program would be as slow as R version. do you have any idea why this problem occur? or did anybody ever have the same problem or can give me some tips to fix it? ps: all datasets are generated through the sampe process(with a single simulation code)
I'm attaching a cleanup of the code using inline -- you can now just
source this in R. I'm not sure, but I think you might be working
under some mistaken assumptions. Hopefully the attached code is more
conceptually clear.
There were a few offhand things that I modified, e.g. maxItr wasn't
used in the attached code. There are a number of things that I
*didn't* change for which there seem to be faster/better methods, but
I'm not sure of your purpose. For example, can any of the in-loop
assignments (qnum, qlength) be moved to the top?
Lastly, are you *just* using rmultinom and Cwhich to roll a K or P
sided die? If so, what about something like as<int>(runif(1, 0, K))?
Even if you're rolling a weighted die, it seems like there's a better
way...
Hopefully this sheds at least a bit of light on your question.
Apologies if I've misunderstood.
-Christian
## Compile/run from R using the following (assuming inline package is
installed):
source('gibbs.R.txt')
fun1(10)
A man, a plan, a cat, a ham, a yak, a yam, a hat, a canal ? Panama!
-------------- next part --------------
#include "SamplingCode.h"
inc = '
//using namespace Rcpp;
int Cwhich(IntegerVector binVec){
int vecsize = binVec.size();
for ( int i =0; i < vecsize; i++)
if (binVec(i) == 1)
return i;
std::cout << "soemthing went wrong, no value equal to one was found !" << std::endl;
return -1;
}
'
#RcppExport SEXP doGibbs(SEXP RmaxItr){//doGibbs <- function(Questions, K, P, gibbs.max.iter, Z.Q=NULL, Y.Q = NULL){
src = '
Environment stats("package:stats");
Function rmultinom = stats["rmultinom"];
Rcpp::RNGScope scope;
// Rcpp::List questions(Rquestions);
unsigned int K =50;
unsigned int P = 50;
NumericVector probs1(P, 1.0/P); //added -- instead of rep
NumericVector probs2(K, 1.0/K); //added
unsigned int maxItr = Rcpp::as<int>(RmaxItr);
Rcpp::List ZQ;
Rcpp::List YQ;
unsigned int unum=300;
std::cout<<unum<<std::endl;
std::cout<<"point1"<<std::endl;
for(unsigned itr=0; itr<maxItr; itr++){ //changed
std::cout<<"itr"<<itr<<std::endl;
for ( unsigned int u=0; u < unum; u++) {
Rcpp::List QZQ;
unsigned int qnum= 100;
IntegerVector tmpu(qnum); //added
//YQ.insert(u, rep(0, qnum ));
for( unsigned int q=0; q<qnum; q++){
//as<Rcpp::IntegerVector>(YQ(u))(q) = Cwhich(rmultinom(1,1, probs1 ));
tmpu(q) = Cwhich(rmultinom(1,1, probs1 ));
unsigned int qlength=200;
IntegerVector tmpq(qlength); //added
//QZQ.insert(q, rep(0, qlength ));
for( unsigned int n=0; n<qlength; n++){
tmpq(n)= Cwhich(rmultinom(1,1, probs2 ));
//as<Rcpp::IntegerVector>(QZQ(q))(n)= Cwhich(rmultinom(1,1, probs2 ));
}
QZQ.insert(q, tmpq);
}
YQ.insert(u, tmpu);
ZQ.insert(u, QZQ);
}
}
'
require(inline)
fun1 = cxxfunction(signature(RmaxItr='integer'), body=src, includes=inc, plugin='Rcpp')
#}