14
$\begingroup$

Observe:

TableForm@
 Map[{First[Timing[RandomVariate[BinomialDistribution[10 #, 1/#]]]],
      First[Timing[RandomVariate[BinomialDistribution[10 #, 1/(# + 1)]]]]} &,
  {1500, 3000, 5000, 10000}]

0.          0.733205
0.          1.809612
0.          3.447622
0.          0.

Note how the timing goes through the roof simply lowering the probability a tiny amount.

If appears for lowish n in Binomial[n,p], when n*p<10, things get strange.

9.01 Windows, would appreciate if others could confirm behavior.

N.b - the times above are on my loungebook, so you might need to use this on a real machine to compare - do you see the same drop?:

TableForm@
 Map[{First[Timing[Do[RandomVariate[BinomialDistribution[10 #, 1/#]], {100}]]],
    First[Timing[Do[RandomVariate[BinomialDistribution[10 #, 1/(# + 1)]], {100}]]]} &,
  {1500, 3000, 5000, 10000}]

Update: This appears to be caused by heuristics in method switching (hazily outlined in the documentation), probably between a lookup/alias and other methods. This can be seen from:

Table[{MaxMemoryUsed[RandomVariate[BinomialDistribution[10 x, 1/(x)]]], 
       MaxMemoryUsed[RandomVariate[BinomialDistribution[10 x, 1/(x + 1)]]]}, 
    {x, {9999, 10000}}]

(* {{3320, 12630024}, {3200, 2440}} *)

Where in this case, x>=10K uses an alternate method.

$\endgroup$
  • $\begingroup$ {0.015625, 0.03125}, {0.`, 0.0625}, {0.`, 0.125}, {0.`, 0.`}} Windows 8.1, Mathematica 10.0.0.0 $\endgroup$ – Sektor Aug 26 '15 at 8:28
  • $\begingroup$ MMA 10.2, Ubuntu, I get {0.03, 0.1, 0.17, 0.0001} for the second column $\endgroup$ – dr.blochwave Aug 26 '15 at 8:33
  • $\begingroup$ MMA 10.2, MBP OS 10.10.5, {{0.000292, 0.027283}, {0.000273, 0.058267}, {0.00024, 0.11855}, {0.000227, 0.000081}} $\endgroup$ – Bob Hanlon Aug 26 '15 at 13:47
  • $\begingroup$ MMA 9.0.1 kubuntu 14.04, {{0.03036, 4.441168}, {0.020735, 11.5795}, {0.025354, 23.53534}, {0.019794, 0.011695}} $\endgroup$ – KennyColnago Aug 26 '15 at 14:17
  • $\begingroup$ What if AbsoluteTiming[] was used instead? $\endgroup$ – J. M. is away Aug 26 '15 at 15:58
6
$\begingroup$

RandomVariate for BinomialDistribution[n,p] changes between methods depending on the value of Min[n*{p,1-p}]. What we're seeing here is that one of those methods is poorly optimized.

Because of this thread, we've made some improvements which should improve speed when Min[n*{p,1-p}]<10. These will be in the next release of Mathematica. We'll also investigate moving to a better algorithm in the future.

$\endgroup$
  • $\begingroup$ RandomVariate functionality in general needs a lot of improvement.This is one place where symbolic capabilities are not needed, so if random number generation can be made compilable (at least for few common distributions, e.g., for those in the Rmath library of R), it would be a great help in simulation work. $\endgroup$ – Asim Aug 28 '15 at 20:55
  • $\begingroup$ Good to hear, thanks for answering. $\endgroup$ – ciao Aug 28 '15 at 23:41
  • $\begingroup$ @ciao looks fixed in v11 $\endgroup$ – dr.blochwave Oct 6 '16 at 7:16
  • $\begingroup$ @blochwave - thanks for the ping! $\endgroup$ – ciao Oct 6 '16 at 8:03
5
$\begingroup$

I get on Mathematica 10.2, Ubuntu 14.04

In[10]:= Map[{First[
    Timing[Do[
      RandomVariate[BinomialDistribution[10 #, 1/#]], {100}]]], 
   First[Timing[
     Do[RandomVariate[
       BinomialDistribution[10 #, 1/(# + 1)]], {100}]]]} &, {1500, 
  3000, 5000, 10000}]

Out[10]= {{0.023484, 2.37428}, {0.012502, 6.22335}, {0.013843, 
  12.4218}, {0.012031, 0.005005}} 

but when I use reals for the probability, we see,

In[11]:= Map[{First[
    Timing[Do[
      RandomVariate[BinomialDistribution[10 #, 1.0/#]], {100}]]], 
   First[Timing[
     Do[RandomVariate[
       BinomialDistribution[10 #, 1.0/(# + 1.0)]], {100}]]]} &, {1500,
   3000, 5000, 10000}]

Out[11]= {{0.013413, 0.231714}, {0.007382, 0.335648}, {0.007443, 
  0.389959}, {0.006538, 0.003143}}
$\endgroup$
5
$\begingroup$

Please edit with your results:

MMa 11.0.0, Ubuntu - blochwave

0.021172    0.019485
0.012286    0.035213
0.012411    0.055937
0.012053    0.005847

MMa 10.0.0.0, Windows 8.1 – Sektor

{0.015625, 0.03125}, 
{0.`,      0.0625}, 
{0.`,      0.125}, 
{0.`,      0.`}}

MMa 10.0.0.0 through MinGW & mintty, Windows 8.1 – Sektor

{0.,   0.03125}, 
{0.,   0.0625}, 
{0.,   0.125}, 
{0.,   0.}}

MMA 10.2, Ubuntu 12.04 - blochwave

          {0.03,
           0.1, 
           0.17, 
           0.0001}

MMA 10.2, MBP OS 10.10.5, - Bob Hanlon

{{0.000292, 0.027283}, 
{0.000273,  0.058267}, 
{0.00024,   0.11855}, 
{0.000227,  0.000081}} 

MMA 9.0.1 kubuntu 14.04, - KennyColnago

{{0.03036, 4.441168},
 {0.020735, 11.5795}, 
 {0.025354, 23.53534}, 
 {0.019794, 0.011695}} 

Mma 9.0.1 WinXP, Belisarius

{{0.,       0.328125}, 
 {0.,       0.750000},
 {0.,       1.484375}, 
 {0.,       0.}}

MMA 10.0.2.0 MBAir OSX 10.10.5 - march

0.029516    3.184197
0.024093    8.040635
0.018399    15.686205
0.023200    0.006712

MMA 9.0.1.0, MBPro Retina, 15-inch, Late 2013, OS 10.9.5 - heropup

{{0.012290, 0.040470},
{0.000434, 0.081408},
{0.000289, 0.162845},
{0.000271, 0.000141}}
$\endgroup$
  • $\begingroup$ Well, interesting. Bug? Or boo-boo algorithmic decision (as in it's switching to something I'd call broken with the suspiciously familiar n*p <10 rule) $\endgroup$ – ciao Aug 26 '15 at 19:15
  • $\begingroup$ @ciao I would vote for the second option ... $\endgroup$ – Dr. belisarius Aug 26 '15 at 19:19
  • $\begingroup$ yep, see update to OP... $\endgroup$ – ciao Aug 28 '15 at 0:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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