稳定分布的随机数？答案

【问题标题】：Stable distribution random numbers?稳定分布的随机数？
【发布时间】：2012-03-10 22:49:58
【问题描述】：

如何在 C# 中生成分布稳定的随机数？ Random 类具有均匀分布。上的许多其他代码互联网显示正态分布。但是我们需要稳定的分布表示无限方差，也就是肥尾分布。

原因是为了生成现实的股票价格。在真实在世界范围内，价格发生巨大变化的可能性远大于正态分布。

有人知道转换随机类输出的 C# 代码吗稳定分布？

编辑：嗯。准确的分布不如确保它会随机生成巨大的 sigma（如至少 20 sigma）那么重要。我们想测试一个交易策略在真正的肥尾分布中的弹性，这正是股票市场价格的表现。

由于 cmets，我刚刚阅读了有关 ZipFian 和 Cauchy 的信息。由于我必须选择，让我们使用 Cauchy 分布，但我也会尝试使用 ZipFian 进行比较。

【问题讨论】：

维基百科不同意您对稳定意味着什么的描述：en.wikipedia.org/wiki/Stable_distribution 例如，它说“正态分布是一组稳定分布。”
我认为你需要找到你需要的确切分布的数学定义（正如@David所说，它可能应该有另一个名字），然后实现你自己的函数将均匀分布转换为一个你需要。
最重要的是你可以学习一些基本的数学。 Fisher 变换（一个指标）进行高斯归一化 - 你所做的是相反的。 Bascialyl 采用更大的随机空间 (0-1000) 来生成数字 0-100 并将更大的区域分配给远离 0 的数字，可以这么说。这实际上只是均匀分布的“高斯”，也是我向愿意从事定量研究的程序提出的问题之一；）以确保他们了解统计数据 101；）
问题应该迁移到http://math.stackexchange.com/吗？
@DavidHeffernan：维基百科页面有些混乱。我认为正态分布是稳定分布的一种非常特殊的情况，因为它是唯一没有肥尾的稳定分布。但是，我已经很久很久没有做统计了，所以我在这里可能错了。

标签： c# random distribution

【解决方案1】：

一般来说，方法是：

选择稳定的肥尾分布。比如说柯西分布。
查找所选分布的分位数函数。

对于 Cauchy 分布，应为 p --> peak + scale * tan( pi * (p - 0.5) )。

现在您有了一种将均匀分布的随机数转换为柯西分布的随机数的方法。

有意义吗？见

http://en.wikipedia.org/wiki/Inverse_transform_sampling

了解详情。

警告：我已经很久很久没有做统计了。

更新：

我非常喜欢这个问题，我刚刚写了博客：见

http://ericlippert.com/2012/02/21/generating-random-non-uniform-data/

我的文章探讨了一些有趣的 Zipfian 分布示例：

http://blogs.msdn.com/b/ericlippert/archive/2010/12/07/10100227.aspx

【讨论】：

这适用于大多数发行版。现在，至于正态分布，那就完全是另一回事了！

【解决方案2】：

如果您对使用 Zipfian 分布（在对科学或社会领域的过程进行建模时经常使用）感兴趣，您可以执行以下操作：

为分布选择您的 k（偏斜）
预计算累积分布的域（这只是一种优化）
通过查找域中最接近的值来为分布生成随机值

示例代码：

List<int> domain = Enumerable.Range(0,1000);  // generate your domain
double skew  = 0.37; // select a skew appropriate to your domain
double sigma = domain.Aggregate(0.0d, (z,x) => x + 1.0 / Math.Pow(z+1, skew));
List<double> cummDist = domain.Select( 
      x => domain.Aggregate(0.0d, (z,y) => z + 1.0/Math.Pow(y, skew) * sigma));

现在您可以通过从域中选择最接近的值来生成随机值：

Random rand = new Random();
double seek = rand.NextDouble();
int searchIndex = cummDist.BinarySearch(seek);
// return the index of the closest value from the distribution domain
return searchIndex < 0 ? (~searchIndex)-1 : searchIndex-1;

当然，您可以通过从映射并从该域返回值的过程中分解出具体化分布域的逻辑来概括整个过程。

【讨论】：

老兄，我很欣赏你优雅的代码，并挖掘了预先计算值的优化。您知道如何将其转化为时间序列股票价格吗？
似乎柯西更适合股票价格随机数。你不同意吗？我会尝试两者，看看哪个看起来更好。
@Wayne：我让那些更聪明、更有经验的人来建议哪种分布更准确地模拟股票价格波动。有很多关于金融市场建模的文章（维基百科只是冰山一角en.wikipedia.org/wiki/…）。最终，可能只有您能够回答哪个分布最合适的问题，具体取决于您的模拟中最重要的是什么。

【解决方案3】：

我面前有 James Gentle 关于这个主题的 Springer 卷，Random Number Generation and Monte Carlo Methods，由我的统计学家妻子提供。它在第 105 页讨论了稳定的家庭：

稳定的分布家族是一个灵活的家族，通常是重尾分布。该族包括一个参数的一个极值处的正态分布和另一个极值处的柯西分布。 Chambers、Mallows 和 Stuck (1976) 给出了一种生成偏离稳定分布的方法。（注意辅助函数 D2 中的常数中的一些错误，用于评估 (e^x-1)/x。）他们的方法在 IMSL 库中使用。对于对称稳定分布，Devroye (1986) 指出，可以通过利用对称稳定与 Fejer-de la Vallee Poissin 分布的关系来开发更快的方法。 Buckle (1995) 展示了如何根据数据模拟稳定分布的参数。

产生偏离一般稳定分布的偏差是困难的。如果您需要这样做，那么我会推荐一个库，例如 IMSL。我不建议你自己尝试。

但是，如果您正在寻找稳定系列中的特定分布，例如Cauchy，那么您可以使用 Eric 描述的方法，称为probability integral transform。只要你能写出封闭形式的分布函数的倒数，你就可以使用这种方法。

【讨论】：

【解决方案4】：

在给定形状参数alpha 和beta 的情况下，以下C# 代码会生成一个遵循稳定分布的随机数。我在知识共享零下将其发布到公共领域。

public static double StableDist(Random rand, double alpha, double beta){
    if(alpha<=0 || alpha>2 || beta<-1 || beta>1)
        throw new ArgumentException();
    var halfpi=Math.PI*0.5;
    var unif=NextDouble(rand);
    while(unif == 0.0)unif=NextDouble(rand);
    unif=(unif - 0.5) * Math.PI;
    // Cauchy special case
    if(alpha==1 && beta==0)
        return Math.Tan(unif);
    var expo=-Math.Log(1.0 - NextDouble(rand));
    var c=Math.Cos(unif);
    if(alpha == 1){
        var s=Math.Sin(unif);
        return 2.0*((unif*beta+halfpi)*s/c -
            beta * Math.Log(halfpi*expo*c/(
                unif*beta+halfpi)))/Math.PI;
    }
    var z=-Math.Tan(halfpi*alpha)*beta;
    var ug=unif+Math.Atan(-z)/alpha;
    var cpow=Math.Pow(c, -1.0 / alpha);
    return Math.Pow(1.0+z*z, 1.0 / (2*alpha))*
        (Math.Sin(alpha*ug)*cpow)*
        Math.Pow(Math.Cos(unif-alpha*ug)/expo, (1.0-alpha) / alpha);
}

private static double NextDouble(Random rand){
    // The default NextDouble implementation in .NET (see
    // https://github.com/dotnet/corert/blob/master/src/System.Private.CoreLib/shared/System/Random.cs)
    // is very problematic:
    // - It generates a random number 0 or greater and less than 2^31-1 in a 
    // way that very slightly biases 2^31-2.
    // - Then it divides that number by 2^31-1.
    // - The result is a number that uses roughly only 32 bits of pseudorandomness, 
    // even though `double` has 53 bits in its significand.
    // To alleviate some of these problems, this method generates a random 53-bit
    // random number and divides that by 2^53.  Although this doesn't fix the bias
    // mentioned above (for the default System.Random), this bias may be of
    // negligible importance for most purposes not involving security.
    long x=rand.Next(0,1<<30);
    x<<=23;
    x+=rand.Next(0,1<<23);
    return (double)x / (double)(1L<<53);
}

另外，我在a separate article中列出了稳定分发的伪代码。

【讨论】：