1
$\begingroup$

Sometimes distribution is too complicated for built-in PDF command, so it would be useful to have approximatePDF[dist,x,n] which creates $n$th order approximation.

Below is one such example of distribution (from here) and working implementation of approximatePDF for $n=2$, can someone make it work for $n>2$?

(* Linear combination of compnents of n-dimensional Dirichlet and n \
integers decaying at rate s.
See http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.46.4629
 *)
makeDist[n_, s_] := Module[{xs, dist, qs, rv},
   xs = Array[x, n - 1];
   dist = DirichletDistribution[ConstantArray[1/2, n]];
   qs = Table[(1 - 1/i)^(2 s) 1/i, {i, 2, n}];
   rv = qs.xs;
   TransformedDistribution[rv, xs \[Distributed] dist]
   ];
dist = makeDist[8, 2];

(* approximate PDF of order 2 *)

approximatePDF[dist_, x_, 2] := Module[{},
   x1 = Cumulant[dist, 1];
   x2 = Cumulant[dist, 2];
   PDF[NormalDistribution[x1, Sqrt[x2]], x]
   ];

pdf = approximatePDF[dist, x, 
  2]; (* PDF[dist,x] is too slow *)
density = 
 Plot[pdf, {x, 0, .1}, PlotStyle -> {Dashed, Gray}, 
  PlotLegends -> {"2nd order"}];
histogram = 
  SmoothHistogram[RandomVariate[dist, 100000], 
   PlotLegends -> {"histogram"}];
Show[density, histogram, PlotRange -> All]

enter image description here

$\endgroup$
2

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.