{"word_count": 122, "character_count": 1014}
Yesterday, I was asked how to write a code to generate a compound Poisson variables, i.e. a series of random variables  where  is a counting random variable (here Poisson disributed) and where the 's are i.i.d (and independent of ), with the convention  when . I came up with the following algorithm, but I was wondering if it was possible to get a better one...&#62;  rcpd=function(n,rN,rX){+  N=rN(n)+  X=rX(sum(N))+  I=as.factor(rep(1:n,N))+  S=tapply(X,I,sum)+  V=as.numeric(S[as.character(1:n)])+  V[is.na(V)]=0+  return(V)}Here, consider - to illustrate - the case where  and ,&#62;  rN.P=function(n) rpois(n,5)&#62;  rX.E=function(n) rexp(n,2)We can generate a sample&#62;  S=rcpd(1000,rN=rN.P,rX=rX.E)and check (using simulation) than &#62; mean(S)[1] 2.547033&#62; mean(rN.P(1000))*mean(rX.E(1000))[1] 2.548309and that &#62; var(S)[1] 2.60393&#62; mean(rN.P(1000))*var(rX.E(1000))++ mean(rX.E(1000))^2*var(rN.P(1000))[1] 2.621376If anyone might think of a faster algorithm, I'd be glad to hear about it...