# Compound Poisson and vectorized computations

[This article was first published on

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

**Freakonometrics - Tag - R-english**, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

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…

> 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 ,

> rN.P=function(n) rpois(n,5) > rX.E=function(n) rexp(n,2)

We can generate a sample

> S=rcpd(1000,rN=rN.P,rX=rX.E)

and check (using simulation) than

> mean(S) [1] 2.547033 > mean(rN.P(1000))*mean(rX.E(1000)) [1] 2.548309

and that

> var(S) [1] 2.60393 > mean(rN.P(1000))*var(rX.E(1000))+ + mean(rX.E(1000))^2*var(rN.P(1000)) [1] 2.621376

If anyone might think of a faster algorithm, I’d be glad to hear about it…

To

**leave a comment**for the author, please follow the link and comment on their blog:**Freakonometrics - Tag - R-english**.R-bloggers.com offers

**daily e-mail updates**about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.