8
$\begingroup$

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?

$\endgroup$
6
  • 1
    $\begingroup$ RandomInteger[PoissonDistribution[#]] & /@ list is slightly faster $\endgroup$ Oct 6, 2015 at 17:27
  • $\begingroup$ 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) $\endgroup$ Oct 6, 2015 at 17:41
  • $\begingroup$ And in fact, another Poisson-based possible duplicate (35433) $\endgroup$ Oct 6, 2015 at 17:44
  • 1
    $\begingroup$ 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. $\endgroup$
    – Ian Hincks
    Oct 6, 2015 at 18:03
  • $\begingroup$ 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 $\endgroup$ Oct 6, 2015 at 18:06

1 Answer 1

5
$\begingroup$

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

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

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.