如何在给定scipy/numpy的拟合曲线情况下最大化直方图的泊松似然?
我在使用python/numpy/scipy的环境中有一些数据,需要把这些数据拟合到一个概率密度函数上。有一种方法是先把数据做成直方图,然后再给这个直方图拟合一条曲线。这里用到的方法是scipy.optimize.leastsq
,它通过最小化(y - f(x))**2
的总和来实现,其中(x,y)在这里指的是直方图的每个小区间的中心和对应的数量。
从统计学的角度来看,这种最小二乘法是为了最大化通过从一个围绕拟合函数的高斯分布中抽样每个小区间的计数而得到这个直方图的可能性。你可以很容易理解这一点:每一项(y-f(x))**2
实际上是-log(gauss(y|mean=f(x)))
,而这些项的总和就是所有小区间的高斯可能性相乘后的对数。
不过,这种方法并不总是准确的:对于我正在研究的统计数据,每个小区间的计数实际上是一个泊松过程的结果,所以我想要最小化所有小区间(x,y)的poisson(y|mean=f(x))
的乘积的对数。泊松分布在f(x)值很大的时候和高斯分布非常接近,但如果我的直方图统计数据不够好,这两者之间的差异就会变得重要,并影响拟合的结果。
1 个回答
0
如果我理解得没错,你有一些数据,想看看某种概率分布是否适合你的数据。
如果是这样的话,你需要用到QQ图。你可以看看这个StackOverflow上的问答。不过那个主要是关于正态分布的,而你需要的是泊松分布的代码。你只需要根据泊松随机函数生成一些随机数据,然后用这些样本进行测试。这里有一个关于泊松分布的QQ图的例子。以下是这个网站上的代码:
#! /usr/bin/env python
from pylab import *
p = poisson(lam=10, size=4000)
m = mean(p)
s = std(p)
n = normal(loc=m, scale=s, size=p.shape)
a = m-4*s
b = m+4*s
figure()
plot(sort(n), sort(p), 'o', color='0.85')
plot([a,b], [a,b], 'k-')
xlim(a,b)
ylim(a,b)
xlabel('Normal Distribution')
ylabel('Poisson Distribution with $\lambda=10$')
grid(True)
savefig('qq.pdf')
show()