       Re: TransformedDistribution

• To: mathgroup at smc.vnet.net
• Subject: [mg116920] Re: TransformedDistribution
• From: Darren Glosemeyer <darreng at wolfram.com>
• Date: Fri, 4 Mar 2011 03:40:58 -0500 (EST)

```On 3/3/2011 4:56 AM, Jagra wrote:
> I have developed some code to generate random variates from the
> product of a LogNormalDistribution and a StableDistribution:
>
> LNStableRV[{\[Alpha]_, \[Beta]_, \[Gamma]_, \[Sigma]_, \[Delta]_},
>
> n_] := Module[{LNRV, SDRV, LNSRV},
>    LNRV = RandomVariate[LogNormalDistribution[Log[\[Gamma]], \[Sigma]],
>       n];
>    SDRV = RandomVariate[
>      StableDistribution[\[Alpha], \[Beta], \[Gamma], \[Sigma]], n];
>    LNRV * SDRV + \[Delta]
>    ]
>
> (* Note the delta serves as a location parameter *)
>
> I think this works fine:
>
> LNStableRV[{1.5, 1, 1, 0.5, 1}, 50000];
> Histogram[%, Automatic, "ProbabilityDensity",
>   PlotRange ->  {{-4, 6}, All}, ImageSize ->  250]
> ListPlot[%%, Joined ->  True, PlotRange ->  All]
>
> Now I'd like to create a TransformedDistribution along the same lines
> so that I can use PDF[], CDF[], etc. on this custom distribution and
> easily do plots and other analysis.
>
> Extrapolating from an example in Documentation Center>
> TransformedDistribution:
>
> \[ScriptCapitalD] =
>    TransformedDistribution[
>     u v, {u \[Distributed] ExponentialDistribution[1/2],
>      v \[Distributed] ExponentialDistribution[1/3]}];
>
> I've tried this:
>
> LogNormalStableDistribution[\[Alpha]_, \[Beta]_, \[Gamma]_, \
> \[Sigma]_, \[Delta]_] := Module[{u, v},
>     TransformedDistribution[
>      u * v + \[Delta], {u \[Distributed]
>        LogNormalDistribution[Log[\[Gamma]], \[Sigma]],
>       v \[Distributed]
>        StableDistribution[\[Alpha], \[Beta], \[Gamma], \[Sigma]]}]
>     ];
>
> \[ScriptCapitalD] = LogNormalStableDistribution[1.5, 1, 1, 0.5, 1]
>
> Which gives me this:
>
> TransformedDistribution[
>   1 + \[FormalX]1 \[FormalX]2, {\[FormalX]1 \[Distributed]
>     LogNormalDistribution[0, 0.5], \[FormalX]2 \[Distributed]
>     StableDistribution[1, 1.5, 1, 1, 0.5]}]
>
> But when I try to plot a PDF of the distribution it never seems to
> finish (granted I haven't let it run more than 10 minutes so far):
>
> Plot[PDF[\[ScriptCapitalD], x], {x, -4, 6}] (* This should plot over
> the same range as the Histogram above *)
>
> So, some questions:
>
> Does my function: LogNormalStableDistribution[] make sense to do this
> kind of thing?
>
> If yes do I:
>
> Just need to let the Plot[] run longer?
> Change it somehow?
> What can I do to make it run faster?
>
> If not:
>
> Do I need to approach this in a different way?
> Use MixtureDistribution?
> Use Something else?
> Any ideas appreciated.
>
> Best,
>
> J
>

Computing a stable PDF value is typically tough to do and requires
numeric integration. Combining a stable with a lognormal makes the
problem more difficult. In this case, it appears that PDF is trying very
hard to compute individual values, but the symbolic methods used by
TransformedDistribution take a long time for each point and are
typically unable to get results.

Random numbers, however, do not require computation of the PDF and are
quite fast. In cases such as this, a good approximation can often be
obtained by generating a bunch of numbers from the distribution and
constructing a distribution from the data. A kernel density estimation
can be one useful method. In the case at hand, you could use
SmoothKernelDistribution to construct such a distribution and use that
to get approximate results for the transform you're after. The default
settings for SmoothKernelDistribution do pretty well to quickly
approximate the distribution from 10^5 data points:

data = RandomVariate[LogNormalStableDistribution[1.5, 1, 1, 0.5, 1], 10^5];

dist = SmoothKernelDistribution[data]

Show[Histogram[data, Automatic, "PDF", PlotRange -> {{-4, 10}, All}],
Plot[PDF[dist, x], {x, -4, 10}, PlotRange -> All]]

EmpiricalDistribution is another possible approximation, but that will
be a discrete distribution rather than a continuous one, and so PDF
values will be 0 at any point not in the data set.

Keep in mind that any approximation method for long-tailed distributions
such as this will tend to have more deviations from the true
distribution the farther you get into the tails of the distribution.
This would be true of smoothing methods, empirical estimates, or methods
based on numeric integration.

Darren Glosemeyer
Wolfram Research

```

• Prev by Date: Re: making something autoexecute before normal execution
• Next by Date: AstronomicalData errors and internet connectivity