用python拟合直方图

Posted 2023-03-12

【中文标题】用python拟合直方图

【英文标题】：Histogram fitting with python

【发布时间】：2016-02-22 00:51:48

【问题描述】：

我一直在冲浪，但没有找到正确的方法来执行以下操作。

我有一个使用 matplotlib 完成的直方图：

hist, bins, patches = plt.hist(distance, bins=100, normed='True')

从图中，我可以看到分布或多或少是指数分布（泊松分布）。考虑到我的 hist 和 bins 数组，我怎样才能做到最佳拟合？

更新

我正在使用以下方法：

x = np.float64(bins) # Had some troubles with data types float128 and float64
hist = np.float64(hist)
myexp=lambda x,l,A:A*np.exp(-l*x)
popt,pcov=opt.curve_fit(myexp,(x[1:]+x[:-1])/2,hist)

但我明白了

---> 41 plt.plot(stats.expon.pdf(np.arange(len(hist)),popt),'-')

ValueError: operands could not be broadcast together with shapes (100,) (2,)

【问题讨论】：

***.com/questions/25828184/… 或这个***.com/questions/6620471/… 的可能重复项？

我不确定您为什么要使用直方图。所有常见的分布都可以通过几个形状/位置参数来固定。这些通常可以从数据本身非常有效地估计出来。

我回答了一个类似的问题yesterday，你甚至可以在那里找到如何使用你自己的拟合模型。您应该将x=(bins[1:]+bins[:-1])/2; y=hist 作为拟合程序的输入。 @cel：对于噪声数据，最小二乘拟合比分布矩的原始估计要可靠得多。(citation needed)

我不太明白你在那里给出了什么。我正在尝试使用param=spy.expon.fit(distance); a = np.linspace(0,bins[-1],1000, dtype='f'); pdf_fitted=spy.expon.pdf(a,loc=param[0],scale=param[1]);plt.plot(a,pdf_fitted,'r-')

不要plt.plot(stats.expon.pdf(np.arange(len(hist)),popt),'-')，而是plt.plot((x[1:]+x[:-1])/2,myexp((x[1:]+x[:-1])/2,*popt),'-')（或者你喜欢的任何x数组）。

【参考方案1】：

您描述的是exponential distribution 的一种形式，并且您想根据数据中观察到的概率密度估计指数分布的参数。而不是使用非线性回归方法（假设残差是高斯分布的），一种正确的方法可以说是 MLE（最大似然估计）。

scipy在其stats库中提供了大量的连续分布，MLE是用.fit()方法实现的。当然指数分布是there：

In [1]:

import numpy as np
import scipy.stats as ss
import matplotlib.pyplot as plt
%matplotlib inline
In [2]:
#generate data 
X = ss.expon.rvs(loc=0.5, scale=1.2, size=1000)

#MLE
P = ss.expon.fit(X)
print P
(0.50046056920696858, 1.1442947648425439)
#not exactly 0.5 and 1.2, due to being a finite sample

In [3]:
#plotting
rX = np.linspace(0,10, 100)
rP = ss.expon.pdf(rX, *P)
#Yup, just unpack P with *P, instead of scale=XX and shape=XX, etc.
In [4]:

#need to plot the normalized histogram with `normed=True`
plt.hist(X, normed=True)
plt.plot(rX, rP)
Out[4]:

您的distance 将在此处替换X。

【讨论】：

Here are all the scipy.stats distributions PDFs with example code.

normed 已弃用，现在称为density