7
$\begingroup$

I've worked for several days trying to get this code to run faster. The imported data is 15MB.

    cyclesperday = 
  CloudGet["https://www.wolframcloud.com/objects/heaneym/\
cyclesperday"];
nc = 10^5;
AbsoluteTiming[
 samplemanycyclesper5years = 
  ParallelTable[Total[RandomChoice[cyclesperday, 1826]], {100*nc}];
 cycledatatofit = Partition[samplemanycyclesper5years, 100];
 samplecycledistributions = 
  ParallelTable[
   EstimatedDistribution[cycledatatofit[[i]], 
    LogNormalDistribution[\[Mu], \[Sigma]]], {i, 1, nc}];
 cyclesamples = 
  Round[ParallelTable[
    RandomVariate[samplecycledistributions[[j]]], {j, 1, nc}]];
 ]

The above code takes about 110 sec to run.

Here is the code revised as suggested below:

cyclesperday = 
  CloudGet["https://www.wolframcloud.com/objects/heaneym/\
cyclesperday"];
nc = 10^5;
rand = Compile[{{cycl, _Real, 1}, {i, _Integer, 0}}, 
   Module[{pos = RandomInteger[{1, Length[cycl]}, i]}, 
    Total[cycl[[pos]]]], RuntimeAttributes -> {Listable}, 
   Parallelization -> True, RuntimeOptions -> "Speed"];
maxLikelihood = 
  Compile[{{values, _Real, 1}}, 
   Module[{\[Mu] = Mean[Log[values]]}, {\[Mu], 
     Sqrt@Mean[(Log[values] - \[Mu])^2]}], 
   RuntimeAttributes -> {Listable}, Parallelization -> True];
AbsoluteTiming[
 samplemanycyclesper5years = rand[cyclesperday, Array[1826 &, 100*nc]];
 cycledatatofit = Partition[samplemanycyclesper5years, 100];
 samplecycledistributions = maxLikelihood[cycledatatofit];
 cyclesamples = 
  Round[ParallelTable[
    RandomVariate[LogNormalDistribution @@ parms], {parms, 
     samplecycledistributions}]];
 ]

The above revised code takes about 107 sec to run.

How can the different components of the code be significantly faster (see below), but the components put together not be?

Improve this question
$\endgroup$

2 Answers 2

Reset to default
11
$\begingroup$

I'm concentrating on the calculation of samplemanycyclesper5years and samplecycledistributions. For the first one, you select 1826 samples randomly and calculate the total. This is done 10^7 times. We can pack the random total into a compiled function that chooses 1826 random integer positions, accesses cyclesperday and calculates the total

rand = Compile[{{cycl, _Real, 1}, {i, _Integer, 0}},
   Module[{pos = RandomInteger[{1, Length[cycl]}, i]},
     Total[cycl[[pos]]]
    ],
   RuntimeAttributes -> {Listable},
   Parallelization -> True,
   RuntimeOptions -> "Speed"
   ];

The parameter i is how many random values of cycl should be totaled. In your case always 1826. Let's test this

rand[cyclesperday, Array[1826 &, 10^5]]; // AbsoluteTiming
(* {0.921493, Null} *)

and compare

ParallelTable[Total[RandomChoice[cyclesperday, 1826]], {10^5}]; // AbsoluteTiming
(* {5.89441, Null} *)

So this needs only 15% of the time your ParallelTable needs. The next step is to do the same for the estimation of the LogNormalDistribution. The estimation of the parameters is actually very simple with a maximum likelihood estimator and you can write this down yourself

maxLikelihood = Compile[{{values, _Real, 1}},
  Module[{μ = Mean[Log[values]]},
   {μ, Sqrt@Mean[(Log[values] - μ)^2]}
   ],
  RuntimeAttributes -> {Listable},
  Parallelization -> True
]

First a quick check:

EstimatedDistribution[cycledatatofit[[10]], 
 LogNormalDistribution[μ, σ]]
(* LogNormalDistribution[7.42205, 0.042639] *)

maxLikelihood[cycledatatofit[[10]]]
(* {7.42205, 0.042639} *)

Excellent. Now let's time it

samplecycledistributions = 
   ParallelTable[EstimatedDistribution[cycledatatofit[[i]], 
     LogNormalDistribution[μ, σ]], {i, 1, nc}]; // AbsoluteTiming
(* {15.8202, Null} *)

and

maxLikelihood[cycledatatofit]; // AbsoluteTiming
(* {0.167895, Null} *)

So this needs only 1% of the original time. Your complete calculation looks like this

nc = 10^5;
samplemanycyclesper5years = 
 rand[cyclesperday, Array[1826 &, 100*nc]];
cycledatatofit = Partition[samplemanycyclesper5years, 100];
samplecycledistributions = maxLikelihood[cycledatatofit];
cyclesamples = 
 Round[ParallelTable[
   RandomVariate[LogNormalDistribution @@ parms], {parms, 
    samplecycledistributions}]];

and I was able to bring it from 654 seconds to 53 seconds. I checked the final histograms and they match perfectly, but please verify each step yourself.

Improve this answer
$\endgroup$
5
  • 1
    $\begingroup$ That is incredibly helpful! How can I learn how to do what you did? $\endgroup$
    – Michael B. Heaney
    Commented Jun 28, 2018 at 2:23
  • 1
    $\begingroup$ It takes some experience because simply throwing ParallelXXX into the code rarely works. Usually, I have a hunch which parts can be improved. In your case, I saw that the things you do in parallel are rather small computations and the highlevel approach with ParallelTable might be improved by only using thread-parallelization and not kernel-parallelization. For the estimation, I knew how Max-Likelihood worked and I guessed that the formulas for your distribution are possibly rather easy and can be compiled. What really helps is when you have at least a clue how you would write functions.. $\endgroup$
    – halirutan
    Commented Jun 28, 2018 at 14:20
  • $\begingroup$ like EstimatedDistribution if you wouldn't have the whole machinery of Mathematica. And then you need to be persistent until you really understood why something is slow and if you can improve it. This takes time that I have already investigated when hitting problems myself, but it is nothing you cannot learn. If you like home-work for the upcoming weekend: Why don't you try to apply the same approach to the last slow RandomVariate part of your code? It is not as easy as my examples, but you would gain a lot of understanding. $\endgroup$
    – halirutan
    Commented Jun 28, 2018 at 14:24
  • $\begingroup$ Sadly, it would probably take me months of work to get good at it, and my work project doesn't allow that time. I'll try to learn it as I go along. $\endgroup$
    – Michael B. Heaney
    Commented Jun 29, 2018 at 16:40
  • $\begingroup$ I am comparing the run times of the two nb's, and they are about the same. Did I make a mistake somewhere? $\endgroup$
    – Michael B. Heaney
    Commented Jul 10, 2018 at 0:27
6
$\begingroup$

Process cyclesperday thus:

cyclesperday = Developer`ToPackedArray[cyclesperday, Real];

Your data gets loaded as a mix of Real and Integer values. To take best advantage of the CPU, the data should all be of the same type and in a packed array. The second argument causes the integer values to be converted to Real.

Then the rest of the code takes about 150 sec. to run. (I cannot tell you how long the original takes on my computer, because it ran out of memory. Packed arrays save memory, too.)

Improve this answer
$\endgroup$
3
  • $\begingroup$ Thanks for your advice, I now use it to store my data. $\endgroup$
    – Michael B. Heaney
    Commented Jun 29, 2018 at 23:12
  • $\begingroup$ If the data were all integers, would it be faster to store it as integers? $\endgroup$
    – Michael B. Heaney
    Commented Jun 29, 2018 at 23:22
  • $\begingroup$ @MichaelB.Heaney On most machines, integers and reals are 64 bits, so it should take the same amount of time to ship the data out to a device (if that's what you mean by "to store it"). If you're converting to some data format, then I'm not really sure which is faster to convert. It may depend on the format chosen. To compute with the data, just to cover another base, at least with the Intel i7, the throughput is nearly the same. Integers may be a little faster, except when dividing(!); and they overflow sooner. $\endgroup$
    – Michael E2
    Commented Jun 29, 2018 at 23:38

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.