1
$\begingroup$
compiled = Compile[{mu, sigma, x}, 
   Module[{}, Log@PDF[NormalDistribution[mu, Abs[sigma]], x]], 
   RuntimeOptions -> "Speed"];

compiled[0.000001, 0.0000444, 1]
->CompiledFunction::cfn: 

it seems that the input values's are too small. does anyone know how to avoid the error in this case?

Added: Replacing the function with reference to comments, this become faster 10x!(in the whole optimization process including other function) Here is the code.

compiled=Compile[{mu, sigma, x}, Module[{}, Log[1 + 1/(E^((x - mu)^2/(2*Abs[sigma]]^2))*(Sqrt[2*Pi]*Abs[sigma]))]], RuntimeOptions -> "Speed"]; 
$\endgroup$
5
$\begingroup$

The problem is that the term Exp[-((x-mu)/sigma)^2/2] in PDF[NormalDistribution[mu, Abs[sigma]], x] will cause an underflow for values of x whose absolute value of x-mu is large compared to sigma; the resulting values are just too small to be representable in machine precision.

compiled2 = Compile[{mu, sigma, x},
   If[Abs[x - mu] < 26. sigma,
    Evaluate[Log@PDF[NormalDistribution[mu, Abs[sigma]], x]],
    0.
    ],
   RuntimeOptions -> "Speed"
   ];

Edit

I overlooked the appearance of Log. One should expect that it will get rid of the Exp in PDF[NormalDistribution[mu, Abs[sigma]], x]. But Mathematica does not know that all variables are real values and that sigma is also positive. The issue can be resolved by telling her explicitly what your assumptions are. Moreover, Blocking the symbols that are needed in the simplification adds at lot of robustness to the code and adding the options RuntimeAttributes -> Listable, Parallelization -> True, allows for fast evaluation on long lists:

compiled3 = Block[{x, mu, sigma},
   With[{code = FullSimplify[
       Log@PDF[NormalDistribution[mu, Abs[sigma]], x],
       {x ∈ Reals, mu ∈ Reals, sigma > 0}
       ]},
    Compile[{{mu, _Real}, {sigma, _Real}, {x, _Real}},
     code,
     RuntimeAttributes -> Listable,
     Parallelization -> True,
     RuntimeOptions -> "Speed"]
    ]
   ];

Here is a usage and timing example:

n = 20000000;
μ = RandomReal[{-0.00001, 0.00001}, n];
σ = RandomReal[{0, 2 0.0000444}, n];
x = RandomReal[{-1, 1}, n];
compiled3[μ, σ, x]; // RepeatedTiming // First

0.990

Also interesting to know: Since all the functions involved are vectorized right from the start, it is actually faster not to compile the code:

Block[{x, mu, sigma},
  f = {mu, sigma, x} \[Function] Evaluate[FullSimplify[
      Log@PDF[NormalDistribution[mu, Abs[sigma]], x],
      {x ∈ Reals, mu ∈ Reals, sigma > 0}
      ]
     ]
  ];
b = f[μ, σ, x]; // RepeatedTiming // First
Max[Abs[(a - b)/a]]

0.869

1.59214*10^-13

By the way:

CompiledFunctionTools`CompilePrint[compiled] showed that there was a call to MainEvaluate in the compiled function. So, you probably would not have experienced any speedup. Evaluate repairs this.

$\endgroup$
  • $\begingroup$ Why Abs[x-mu]<26. sigma? $\endgroup$ – Xminer Dec 16 '18 at 13:02
  • 1
    $\begingroup$ Because Exp[-27.^2] causes an underflow. Moreover, Exp[-26.^2/2] has already magnitude $10^{-147}$, so setting it to 0. will most likely cause no harm. $\endgroup$ – Henrik Schumacher Dec 16 '18 at 14:06
  • 2
    $\begingroup$ Sqrt[-Log[$MinMachineNumber]] evaluates to 26.6157, so 26. is actually a pretty good cutoff value. $\endgroup$ – J. M. will be back soon Mar 28 at 16:27

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.