# Speeding Up RandomVariate

Consider the following code where we first take 1000 independent samples from a Poisson distribution, and then take 1000 independent samples from 1000 different Poisson distributions:

list = RandomVariate[NormalDistribution[0.5, 0.1],1000];
RandomVariate[PoissonDistribution[0.5], 1000] // Timing // First
RandomVariate[PoissonDistribution[#]] & /@ list // Timing // First


which on my machine outputs something like 0.000343 and 0.137541, respectively. Note that it obviously doesn't matter if list contains different numbers or identical numbers.

Clearly RandomVariate has some kind of overhead. Is there a way to improve the efficiency of sampling from a list of distributions that vary only in their parameters?

• RandomInteger[PoissonDistribution[#]] & /@ list is slightly faster – Dr. belisarius Oct 6 '15 at 17:27
• Poisson-specific possible duplicate: (56180) and general-distribution possible duplicate (75303). Note also that attempting to compile with a Poisson distribution won't work, as per (1124) – dr.blochwave Oct 6 '15 at 17:41
• And in fact, another Poisson-based possible duplicate (35433) – dr.blochwave Oct 6 '15 at 17:44
• Thanks. The last one you mentioned seems the most like a duplicate, and there are some good answers there, but I am most interested in a distribution-independent solution, which does not seem to be asked or answered elsewhere. – Ian Hincks Oct 6 '15 at 18:03
• Well (75303) is probably closest - for a distribution-independent solution, you can turn either to a C++ solution, or the clever method presented by @ciao – dr.blochwave Oct 6 '15 at 18:06

## 1 Answer

Per your example case:

(* setup (you'll pay time here, but once done, it's done) *)
dist = ParameterMixtureDistribution[PoissonDistribution[m],
m \[Distributed] NormalDistribution[1/2, 1/10]];

pdf = PDF[dist, z];

(* pick some end-point with negligible probability *)
(* Here, using 30 truncates tail with total p of ~6.34*10^-40 *)
dc = Rule @@ Transpose@Table[{pdf /. z -> zz, zz}, {zz, 0, 30}];

(* use it and compare *)
result = RandomChoice[dc, 1000] // Timing // First

list = RandomVariate[NormalDistribution[0.5, 0.1], 1000];
resultold=RandomVariate[PoissonDistribution[#]] & /@ list // Timing // First


0.

0.842405

Timings on netbook...

• I miss the cigar lounge – Dr. belisarius Oct 7 '15 at 3:25
• @belisariusisforth: there as I type - will take a puff in your honor ;-} – ciao Oct 7 '15 at 5:33
• PDF[dist, z] gives Undefined in 10.3.1. – Kuba Jan 31 '16 at 15:01