4
$\begingroup$

I am obtaining slow gains in compiling a function, and would like to ask if there is something I am obviously missing.

The function I intend to speed up has the form:

Times@@Table[
        PDF[NormalDistribution[.2, 1], r]  / CDF[NormalDistribution[.2, 1],r],
        {r, 1, 100000}]

where in large evaluations the number of elements in the table, and the parameters of the Normal distribution, may change.

This way of evaluating will return a result with a timing of around 8.1 seconds.

I am compiling this function as:

t1 = Compile[{{r, _Real, 0}},
      PDF[NormalDistribution[.2, 1], r]  / CDF[NormalDistribution[.2, 1], r], 
      RuntimeAttributes -> {Listable}];

t2 = Compile[{{rr, _Real, 1}},
      Times@@t1[rr]];

and then invoking

t2[Table[r, {r, 1, 100000}]

In this case, I get a timing of around 6 seconds - an improvement, but not what I expected. I don't see any MainEvaluate calls. Is there anything obvious I am missing?

$\endgroup$
6
  • 1
    $\begingroup$ Wait a minute, you don't see MainEvaluate in CompilePrint[t1]? Did these functions recently become compilable...? I'm on 10.0.2 and CompilePrint[t1] has 5 lines for me, two of which are MainEvaluate :p $\endgroup$ Commented Dec 5, 2016 at 22:03
  • 1
    $\begingroup$ @MariusLadegårdMeyer, same result on 11.0.1 - I think fred is mistaken. $\endgroup$ Commented Dec 5, 2016 at 22:06
  • $\begingroup$ @SimonWoods, thanks for verifying. fred, that explains the lack of speedup of course. $\endgroup$ Commented Dec 5, 2016 at 22:10
  • $\begingroup$ Hi all - you are actually both right, my apologies. If I type "CompiledFunctionTools``CompilePrint[t1]" I only see the instruction; If I call Needs["CompiledFunctionTools"] first, I see them. My inexperience on Compile hits again. [the inverted accent is not showing well in the in-line code above] $\endgroup$
    – Fred
    Commented Dec 5, 2016 at 22:11
  • 1
    $\begingroup$ You should also note that compiled code can only deal with machine numbers and your result is 72 trillion orders of magnitude smaller than the smallest machine real. $\endgroup$ Commented Dec 5, 2016 at 22:26

2 Answers 2

5
$\begingroup$

Since the PDF and CDF of a normal distribution can be expressed analytically as a function of Erfc, you can already speed up your calculation by avoiding the PDF and CDF calculation / lookup at every step:

f[mu_, sigma_] = 
 PDF[NormalDistribution[mu, sigma], r] / CDF[NormalDistribution[mu, sigma], r]

result

Times @@ Table[f[.2, 1.], {r, 1, 100000}]

(*Out: 3.3929765235*10^-72383065132319 *)

This is already significantly faster than the compiled version:

Times @@ Table[f[.2, 1.], {r, 1, 100000}]; // RepeatedTiming

(* Out: {2.7, Null} *)

An alternative formulation, using the Listable properties of the functions involved, would be:

f2[mu_, sigma_][r_] = 
 PDF[NormalDistribution[mu, sigma], r] / CDF[NormalDistribution[mu, sigma], r]

Times @@ f2[0.2, 1][Range[100000]]; // RepeatedTiming

(* Out: {1.164, Null} *)

Of course, the results are all the same:

Times @@ Table[
  PDF[NormalDistribution[.2, 1], r]/CDF[NormalDistribution[.2, 1], r],
  {r, 1, 100000}
  ];
Times @@ f2[0.2, 1][Range[100000]];
Times @@ Table[f[.2, 1.], {r, 1, 100000}];

%%% == %% == %

(* Out: True *)
$\endgroup$
2
  • $\begingroup$ This - plus the comments above on the MainEvaluate calls - is very helpful. I have embedded everything into a single compiled function that makes is very fast now. Will try to update the code above. $\endgroup$
    – Fred
    Commented Dec 5, 2016 at 22:35
  • $\begingroup$ @fred Something probably got lost in your code when copy-pasting. Try to enclose the code in comments inside backticks like this. I would also urge you to write your own answer combining the insight you gained. $\endgroup$
    – MarcoB
    Commented Dec 5, 2016 at 22:41
5
$\begingroup$

Following the suggestions above, here's my improved version.

I have found that coding Erfc directly eliminates the MainEvaluate and makes everything remarkably fast. This my improved code:

t3 = Compile[{{rr, _Real, 1}},
        Times @@ Table[
              (Exp[-((.2 + rr[[f]])/1)^2] ((2/\[Pi])^(1/2))  )/
                      (1 Erfc[(.2 - rr[[f]])/(1 2^0.5)]),
          {f, 1, 100000}]
  ];
$\endgroup$
2
  • $\begingroup$ Fred, Thank you for sharing your solution. Unfortunately, however, I can't seem to get your code to work. Can you show a use example for t3? $\endgroup$
    – MarcoB
    Commented Dec 5, 2016 at 23:34
  • $\begingroup$ Sure, I should have been more clear. In my context, the tensor rr has the same length as the table inside the Compile function. As an example, the instructions x=Table[0.2,{i,1,100000}] and t3[x] work. $\endgroup$
    – Fred
    Commented Dec 6, 2016 at 2:49

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.