【问题标题】:Fitting a gamma distribution with (python) Scipy用(python)Scipy拟合伽马分布
【发布时间】:2011-02-23 04:22:01
【问题描述】:

谁能帮我在 python 中拟合伽马分布?好吧,我有一些数据:X 和 Y 坐标,我想找到适合这个分布的 gamma 参数......在Scipy doc 中,事实证明确实存在 fit 方法但我不知道如何使用它:s.. 首先,参数“数据”必须采用哪种格式,以及我如何提供第二个参数(参数),因为这就是我正在寻找的?

【问题讨论】:

    标签: python scipy distribution gamma-distribution


    【解决方案1】:

    生成一些伽马数据:

    import scipy.stats as stats    
    alpha = 5
    loc = 100.5
    beta = 22
    data = stats.gamma.rvs(alpha, loc=loc, scale=beta, size=10000)    
    print(data)
    # [ 202.36035683  297.23906376  249.53831795 ...,  271.85204096  180.75026301
    #   364.60240242]
    

    这里我们将数据拟合到伽马分布:

    fit_alpha, fit_loc, fit_beta=stats.gamma.fit(data)
    print(fit_alpha, fit_loc, fit_beta)
    # (5.0833692504230008, 100.08697963283467, 21.739518937816108)
    
    print(alpha, loc, beta)
    # (5, 100.5, 22)
    

    【讨论】:

    • 非常感谢!但是为什么一开始就创建了变量 x 呢?
    • 啊,看来我的消息来得太晚了。再次非常感谢你;)
    • scipy.stats 使用最大似然估计进行拟合,因此您需要传递原始数据而不是 pdf/pmf (x, y)
    • 请注意,beta 用于表示分布的“率”参数,它是“形状”的倒数。
    【解决方案2】:

    我对 ss.gamma.rvs 函数不满意,因为它可以生成负数,而 gamma 分布应该没有。所以我通过期望值 = mean(data) 和 variance = var(data) 拟合样本(详见*)并编写了一个函数,可以在没有 scipy 的情况下产生伽马分布的随机样本(我发现很难正确安装,旁注):

    import random
    import numpy
    
    data = [6176, 11046, 670, 6146, 7945, 6864, 767, 7623, 7212, 9040, 3213, 6302, 10044, 10195, 9386, 7230, 4602, 6282, 8619, 7903, 6318, 13294, 6990, 5515, 9157]
    
    # Fit gamma distribution through mean and average
    mean_of_distribution = numpy.mean(data)
    variance_of_distribution = numpy.var(data)
    
    def gamma_random_sample(mean, variance, size):
        """Yields a list of random numbers following a gamma distribution defined by mean and variance"""
        g_alpha = mean*mean/variance
        g_beta = mean/variance
        for i in range(size):
            yield random.gammavariate(g_alpha,1/g_beta)
    
    # force integer values to get integer sample
    grs = [int(i) for i in gamma_random_sample(mean_of_distribution,variance_of_distribution,len(data))]
    
    print("Original data: ", sorted(data))
    print("Random sample: ", sorted(grs))
    
    # Original data: [670, 767, 3213, 4602, 5515, 6146, 6176, 6282, 6302, 6318, 6864, 6990, 7212, 7230, 7623, 7903, 7945, 8619, 9040, 9157, 9386, 10044, 10195, 11046, 13294]
    # Random sample:  [1646, 2237, 3178, 3227, 3649, 4049, 4171, 5071, 5118, 5139, 5456, 6139, 6468, 6726, 6944, 7050, 7135, 7588, 7597, 7971, 10269, 10563, 12283, 12339, 13066]
    

    【讨论】:

    • “参见*了解详细信息”非常通用。您应该添加了一个特定的链接。
    【解决方案3】:

    如果您想要一个长示例,包括有关估计或修复分发支持的讨论,那么您可以在https://github.com/scipy/scipy/issues/1359 和链接的邮件列表消息中找到它。

    在 fit 期间修复参数(例如位置)的初步支持已添加到 scipy 的主干版本中。

    【讨论】:

      【解决方案4】:

      OpenTURNS 使用GammaFactory 类有一个简单的方法。

      首先,让我们生成一个样本:

      import openturns as ot
      gammaDistribution = ot.Gamma()
      sample = gammaDistribution.getSample(100)
      

      然后为它拟合一个 Gamma:

      distribution = ot.GammaFactory().build(sample)
      

      那么我们就可以绘制Gamma的PDF了:

      import openturns.viewer as otv
      otv.View(distribution.drawPDF())
      

      产生:

      有关此主题的更多详细信息,请访问:http://openturns.github.io/openturns/latest/user_manual/_generated/openturns.GammaFactory.html

      【讨论】:

        【解决方案5】:

        1):“data”变量可以是python列表或元组的格式,也可以是numpy.ndarray,可以通过以下方式获取:

        data=numpy.array(data)
        

        上述行中的第二个数据应该是一个列表或元组,包含您的数据。

        2:“参数”变量是第一个猜测,您可以选择将其提供给拟合函数作为拟合过程的起点,因此可以省略。

        3:关于@mondano 答案的注释。使用矩(均值和方差)来计算 gamma 参数对于较大的形状参数(alpha>10)是相当好的,但对于较小的 alpha 值可能会产生较差的结果(请参阅大气科学中的统计方法 Wilks 和 THOM,HCS,1958 年:关于伽马分布的注释。Mon. Wea. Rev., 86, 117–122。

        在这种情况下,使用 scipy 模块中实现的最大似然估计器被认为是更好的选择。

        【讨论】: