4
$\begingroup$

I have the defined a custom distribution, which is a Gaussian where the exponent can be any real n greater than 0 (set to actually be >0.5 in use), i.e.

ProbabilityDistribution[
  (2^(-((2^n*(Abs[x])^n)/FWHM^n))), {x, -Infinity, Infinity}, 
  Method -> "Normalize"
]

Under this normalisation, the part of the Piecewise function I need simplifies down to (when FullSimplify is used):

CustomFunction[FWHM_, n_] = 
  ProbabilityDistribution[
    (2^(-1 - 2^n FWHM^-n Abs[x]^n) (2^n FWHM^-n)^(1/n) Log[2]^(1/n))/Gamma[1 + 1/n](*Re[n]>0*),
    {x, -Infinity, Infinity}
  ]

I wish to generate N values from this distribution, which I normally would do using RandomVariate:

RandomVariate[CustomFunction[FWHM, n] /. {FWHM -> 100, n -> 0.5}, 100(*N*)]

However I have noticed that for some values for FWHM and n this is very slow. For example:

RandomVariate[
  CustomFunction[FWHM, n] /. {FWHM -> 5235.315151`, n -> 0.534124`},
  1000
]; // AbsoluteTiming

takes >70s to run on my PC, whereas other numbers take less than 1s.

Is there any reason for this, or a way to speed up the RandomVariate?

$\endgroup$
1
  • 1
    $\begingroup$ I can't reproduce the behavior you see in MMA 12.3.1 on Win10-64. Your last call to RandomVariate takes 0.07s on my laptop, which is not particularly performant. It takes roughly 2s to get one million variates. Maybe clear everything and restart? $\endgroup$
    – MarcoB
    Sep 9, 2022 at 15:17

1 Answer 1

4
$\begingroup$
$Version

(* "13.1.0 for Mac OS X x86 (64-bit) (June 16, 2022)" *)

Clear["Global`*"]

dist[FWHM_, n_] =
  ProbabilityDistribution[
   (2^(-((2^n*(Abs[x])^n)/FWHM^n))),
   {x, -Infinity, Infinity},
   Method -> "Normalize",
   Assumptions -> {n > 0, FWHM > 0}];

PDF[dist[FWHM, n], x]

(* (2^(-2^n FWHM^-n Abs[x]^n) Log[2]^(1/n))/(FWHM Gamma[1 + 1/n]) *)

CustomFunction[FWHM_?Positive, n_?Positive] := 
    Evaluate@ProbabilityDistribution[
       PDF[dist[FWHM, n], x], 
       {x, -Infinity, Infinity}];

SeedRandom[1234];

(data1 = RandomVariate[CustomFunction[FWHM, n] /. 
      {FWHM -> 100, n -> 0.5}, 1000]); //
 AbsoluteTiming

(* {0.060385, Null} *)

Show[Histogram[data1, Automatic, "PDF"],
 Plot[PDF[CustomFunction[FWHM, n] /.
    {FWHM -> 100, n -> 0.5}, x],
  {x, Min@data1, Max@data1},
  PlotRange -> All]]

enter image description here

Sometimes the data range can be excessive (although the timing is comparable).

SeedRandom[1234];

(data2 = RandomVariate[CustomFunction[FWHM, n] /. 
      {FWHM -> 5235.315151`, n -> 0.534124`}, 1000]); // 
 AbsoluteTiming

(* {0.063217, Null} *)

MinMax@data2

(* {-271327., 4.15334*10^8} *)

Using a different seed

SeedRandom[123];

(data3 = RandomVariate[CustomFunction[FWHM, n] /. 
      {FWHM -> 5235.315151`, n -> 0.534124`}, 1000]); // 
 AbsoluteTiming

(* {0.061821, Null} *)

MinMax@data3

(* {-245766., 356148.} *)

Show[Histogram[data3, Automatic, "PDF"],
 Plot[PDF[CustomFunction[FWHM, n] /.
    {FWHM -> 5235.315151`, n -> 0.534124`}, x],
  {x, Min@data3, Max@data3},
  PlotRange -> All]]

enter image description here

$\endgroup$
2
  • $\begingroup$ Thanks for this, seems to speed it up. Is there a reason why embedding within Evaluate@ProbabilityDistribution[PDF[dist[FWHM, n], x],] helps speed things up? $\endgroup$
    – Xyive
    Sep 12, 2022 at 12:23
  • 1
    $\begingroup$ The revised CustomFunction restricts its arguments to being numeric (Positive) so it never tries to evaluate symbolically. Because of the argument restrictions, it must be defined with SetDelayed rather than Set and the Evaluate is necessary to do the one-time evaluation of PDF[dist[FWHM, n], x]. $\endgroup$
    – Bob Hanlon
    Sep 12, 2022 at 12:41

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.