MathGroup Archive 2011

[Date Index] [Thread Index] [Author Index]

Search the Archive

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
  • Previous by thread: TransformedDistribution
  • Next by thread: Is a conditional evaluation of sections of a notebook possible?