Mathematica 中的变换分布

Question

我开发了一些代码来从 LogNormalDistribution 和 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 *)

我认为这很好用：

LNStableRV[{1.5, 1, 1, 0.5, 1}, 50000];
Histogram[%, Automatic, "ProbabilityDensity",
          PlotRange -> {{-4, 6}, All}, ImageSize -> 250]
ListPlot[%%, Joined -> True, PlotRange -> All]

现在我想创建一个 TransformedDistribution 沿着同样的路线，这样我就可以在这个自定义分布上使用 PDF[]、CDF[] 等，并轻松地进行绘图和其他分析。

从文档中心中的示例推断→TransformedDistribution：

\[ScriptCapitalD] =
  TransformedDistribution[
   u v, {u \[Distributed] ExponentialDistribution[1/2],
    v \[Distributed] ExponentialDistribution[1/3]}];

我尝试过这个：

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]

这给了我这个：

TransformedDistribution[
 1 + \[FormalX]1 \[FormalX]2, {\[FormalX]1 \[Distributed]
   LogNormalDistribution[0, 0.5], \[FormalX]2 \[Distributed]
   StableDistribution[1, 1.5, 1, 1, 0.5]}]

但是当我尝试绘制分布的PDF时，它似乎永远不会完成（当然我没有让它运行超过一分钟或两分钟）：

Plot[PDF[\[ScriptCapitalD], x], {x, -4, 6}] (* This should plot over the same range as the Histogram above *)

所以，有一些问题：

我的函数 LogNormalStableDistribution[] 做这种事情有意义吗？

如果是的话，我：

• 只需要让 Plot[] 运行 更长？
• 以某种方式改变它？
• 我该怎么做才能做到 跑得更快？

如果不是：

• 我需要以不同的方式处理这个问题吗？
• 使用混合分布？
• 使用其他东西？

Answer 1

您使用转换分布的方法很好，但由于分布的 PDF 不以封闭形式存在，因此 PDF[TransformedDistribution[..],x] 不是可行的方法，对于每个 x 将应用符号求解器。最好对你的发行版进行修改以得到 pdf 格式。令 X 为 LogNormal-Stable 随机变量。那么

CDF[LogNormalStableDistribution[params], x] == Probability[X <= x] 

But X==U*V + delta 因此 X<=x 转换为 V<=(x-delta)/U。这给出了

LogNormalStableCDF[{\[Alpha]_, \[Beta]_, \[Gamma]_, \[Sigma]_, \
\[Delta]_}, x_Real] := 
 Block[{u}, 
  NExpectation[
   CDF[StableDistribution[\[Alpha], \[Beta], \[Gamma], \[Sigma]], (x \
- \[Delta])/u], 
   u \[Distributed] LogNormalDistribution[Log[\[Gamma]], \[Sigma]]]]

关于 x 的微分，我们得到 PDF：

LogNormalStablePDF[{\[Alpha]_, \[Beta]_, \[Gamma]_, \[Sigma]_, \
\[Delta]_}, x_Real] := 
 Block[{u}, 
  NExpectation[
   PDF[StableDistribution[\[Alpha], \[Beta], \[Gamma], \[Sigma]], (x \
- \[Delta])/u]/u, 
   u \[Distributed] LogNormalDistribution[Log[\[Gamma]], \[Sigma]]]]

使用它，这里是图