【问题标题】:Fit a distribution to a Counter in scipy将分布拟合到 scipy 中的计数器
【发布时间】:2014-03-19 04:12:21
【问题描述】:

我有一个 collections.Counter 对象,其中包含不同值的出现次数,如下所示:

1:193260
2:51794
3:19112
4:9250
5:6486

如何在 scipy 中为这些数据拟合概率分布? scipy.stats.expon.fit() 似乎想要一个数字列表。用 193260 [1]s、51794 [2]s 等创建一个列表似乎很浪费。有没有更优雅或更有效的方法?

【问题讨论】:

标签: python scipy distribution numerical-methods


【解决方案1】:

看起来 scipy.stats.expon.fit 基本上是 scipy.optimize.minimize 的一个小包装器,它首先创建一个函数来计算 neg-log-likelihood,然后使用 scipy.optimize.minimize 来适应pdf参数。

所以,我认为您需要在这里编写自己的函数来计算计数器对象的负对数似然,然后自己调用 scipy.optimize.minimize。

更具体地说,scipy 在这里定义了 expon 'scale' 参数 http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.expon.html

所以,pdf 是:

pdf(x) = 1 / scale * exp ( - x / scale)

所以,取两边的对数我们得到:

log_pdf(x) = - log(scale) - x / scale

因此,您的 counter 对象中所有内容的负对数相似性将是:

def neg_log_likelihood(scale):
    total = 0.0
    for x, count in counter.iteritems():
       total += (math.log(scale) + x / scale) * count
    return total

这是一个尝试这个的程序。

import scipy.stats
import scipy.optimize
import math
import collections

def fit1(counter):
    def neg_log_likelihood(scale):
        total = 0.0
        for x, count in counter.iteritems():
           total += (math.log(scale) + x / scale) * count
        return total

    optimize_result = scipy.optimize.minimize(neg_log_likelihood, [1.0])
    if not optimize_result.success:
        raise Exception(optimize_result.message)
    return optimize_result.x[0]

def fit2(counter):
    data = []
    # Create an array where each key is repeated as many times
    # as the value of the counter.
    for x, count in counter.iteritems():
        data += [x] * count
    fit_result = scipy.stats.expon.fit(data, floc = 0)
    return fit_result[-1]    

def test(): 
    c = collections.Counter()
    c[1] = 193260
    c[2] = 51794
    c[3] = 19112
    c[4] = 9250
    c[5] = 6486

    print "fit1 'scale' is %f " % fit1(c)
    print "fit2 'scale' is %f " % fit2(c)

test()

这是输出:

fit1 'scale' is 1.513437 
fit2 'scale' is 1.513438 

【讨论】:

    猜你喜欢
    • 2011-10-01
    • 2018-02-08
    • 2016-01-14
    • 1970-01-01
    • 2013-07-03
    • 2021-05-13
    • 2018-08-21
    • 1970-01-01
    • 2011-02-23
    相关资源
    最近更新 更多