Подгонка гистограммы с питоном - PullRequest
36 голосов
/ 18 октября 2011

У меня есть гистограмма

H=hist(my_data,bins=my_bin,histtype='step',color='r')

Я вижу, что форма почти гауссова, но я бы хотел добавить эту гистограмму к гауссовой функции и вывести значение среднего и сигмы, которые я получаю. Вы можете мне помочь?

Ответы [ 4 ]

61 голосов
/ 24 декабря 2011

Вот пример работы с py2.6 и py3.2:

from scipy.stats import norm
import matplotlib.mlab as mlab
import matplotlib.pyplot as plt

# read data from a text file. One number per line
arch = "test/Log(2)_ACRatio.txt"
datos = []
for item in open(arch,'r'):
    item = item.strip()
    if item != '':
        try:
            datos.append(float(item))
        except ValueError:
            pass

# best fit of data
(mu, sigma) = norm.fit(datos)

# the histogram of the data
n, bins, patches = plt.hist(datos, 60, normed=1, facecolor='green', alpha=0.75)

# add a 'best fit' line
y = mlab.normpdf( bins, mu, sigma)
l = plt.plot(bins, y, 'r--', linewidth=2)

#plot
plt.xlabel('Smarts')
plt.ylabel('Probability')
plt.title(r'$\mathrm{Histogram\ of\ IQ:}\ \mu=%.3f,\ \sigma=%.3f$' %(mu, sigma))
plt.grid(True)

plt.show()

enter image description here

25 голосов
/ 20 марта 2013

Вот пример, который использует scipy.optimize для размещения нелинейных функций, таких как гауссиана, даже когда данные находятся в гистограмме, которая не имеет хорошего ранжирования, так что простая средняя оценка не удалась Константа смещения также приведет к сбою простой нормальной статистики (просто удалите p [3] и c [3] для простых гауссовских данных).

from pylab import *
from numpy import loadtxt
from scipy.optimize import leastsq

fitfunc  = lambda p, x: p[0]*exp(-0.5*((x-p[1])/p[2])**2)+p[3]
errfunc  = lambda p, x, y: (y - fitfunc(p, x))

filename = "gaussdata.csv"
data     = loadtxt(filename,skiprows=1,delimiter=',')
xdata    = data[:,0]
ydata    = data[:,1]

init  = [1.0, 0.5, 0.5, 0.5]

out   = leastsq( errfunc, init, args=(xdata, ydata))
c = out[0]

print "A exp[-0.5((x-mu)/sigma)^2] + k "
print "Parent Coefficients:"
print "1.000, 0.200, 0.300, 0.625"
print "Fit Coefficients:"
print c[0],c[1],abs(c[2]),c[3]

plot(xdata, fitfunc(c, xdata))
plot(xdata, ydata)

title(r'$A = %.3f\  \mu = %.3f\  \sigma = %.3f\ k = %.3f $' %(c[0],c[1],abs(c[2]),c[3]));

show()

Выход:

A exp[-0.5((x-mu)/sigma)^2] + k 
Parent Coefficients:
1.000, 0.200, 0.300, 0.625
Fit Coefficients:
0.961231625289 0.197254597618 0.293989275502 0.65370344131

gaussian plot with fit

1 голос
/ 01 марта 2019

Начиная с Python 3.8, стандартная библиотека предоставляет объект NormalDist в составе модуля statistics.

Объект NormalDist можетбыть построенным из набора данных с помощью метода NormalDist.from_samples и обеспечивает доступ к его среднее (NormalDist.mean) и стандартное отклонение (NormalDist.stdev):

from statistics import NormalDist

# data = [0.7237248252340628, 0.6402731706462489, -1.0616113628912391, -1.7796451823371144, -0.1475852030122049, 0.5617952240065559, -0.6371760932160501, -0.7257277223562687, 1.699633029946764, 0.2155375969350495, -0.33371076371293323, 0.1905125348631894, -0.8175477853425216, -1.7549449090704003, -0.512427115804309, 0.9720486316086447, 0.6248742504909869, 0.7450655841312533, -0.1451632129830228, -1.0252663611514108]
norm = NormalDist.from_samples(data)
# NormalDist(mu=-0.12836704320073597, sigma=0.9240861018557649)
norm.mean
# -0.12836704320073597
norm.stdev
# 0.9240861018557649
1 голос
/ 28 сентября 2016

Вот еще одно решение, использующее только пакеты matplotlib.pyplot и numpy.Работает только для гауссовой примерки.Он основан на оценке максимального правдоподобия и уже упоминался в этой теме .Вот соответствующий код:

# Python version : 2.7.9
from __future__ import division
import numpy as np
from matplotlib import pyplot as plt

# For the explanation, I simulate the data :
N=1000
data = np.random.randn(N)
# But in reality, you would read data from file, for example with :
#data = np.loadtxt("data.txt")

# Empirical average and variance are computed
avg = np.mean(data)
var = np.var(data)
# From that, we know the shape of the fitted Gaussian.
pdf_x = np.linspace(np.min(data),np.max(data),100)
pdf_y = 1.0/np.sqrt(2*np.pi*var)*np.exp(-0.5*(pdf_x-avg)**2/var)

# Then we plot :
plt.figure()
plt.hist(data,30,normed=True)
plt.plot(pdf_x,pdf_y,'k--')
plt.legend(("Fit","Data"),"best")
plt.show()

и здесь - это выход.

...