Ein besseres Gaussian zu den Datenpunkten anbringen?

Question

Ich versuche, einen Gaussian an eine Reihe von Datenpunkten anzupassen, die einer Gaußverteilung zu folgen scheinen. Ich habe bereits eine Menge Möglichkeiten dafür überprüft, aber die meisten von ihnen verstehe ich nicht wirklich. Ich fand jedoch eine Lösung, die zu funktionieren scheint, aber die tatsächliche Anpassung, die ich erhalte, sieht nicht viel mehr nach einem Gaussian aus als meine Datenpunkte.

Hier ist mein Code:

import numpy as np 
import matplotlib.pyplot as plt 
from scipy import asarray as ar, exp, sqrt 
from scipy.optimize import curve_fit 


angles = [-8, -6, -4, -2, 0, 2, 4, 6, 8] 
data = [99, 610, 1271, 1804, 1823, 1346, 635, 125, 24] 
angles = ar(angles) 
data = ar(data) 

n = len(x) 
mean = sum(data*angles)/n 
sigma = sqrt(sum(data*(angles-mean)**2)/n) 

def gaus(x,a,mu,sigma): 
    return a*exp(-(x-mu)**2/(2*sigma**2)) 

popt,pcov = curve_fit(gaus,angles,data,p0=[0.18,mean,sigma]) 


fig = plt.figure() 
plt.plot(angles, data, "ob", label = "Measured") 
plt.plot(angles,gaus(angles,*popt),'r',label='Fit') 
plt.xlim(-10, 10) 
plt.ylim(0, 2000) 
plt.xticks(angles) 
plt.title("$^{137}$Cs Zero Point") 
plt.xlabel("Angle [$^\circ$]") 
plt.ylabel("662 keV-Photon Count") 
plt.grid() 
plt.legend() 
plt.show() 

Dies ist die Ausgabe, die es erzeugt:

Wie Sie sehen können, ist der Sitz eine schöne und symmetrische „echte“ nicht beschreiben Gaußsche . Gibt es eine Möglichkeit, einen "besseren" Gaussian zu bekommen oder ist das so gut wie es geht?

Vielen Dank!

Answer 1

Der beste Weg ist, einfach den Mittelwert und die Varianz der Punkte zu verwenden. Ich meine, wenn Sie Zugriff auf die zugrunde liegenden Daten haben, die dieses Histogramm erzeugt haben, dann sollten Sie den Mittelwert und die Varianz unter Verwendung der Funktionen mean und var berechnen.

Das Histogramm ist nur eine visuelle Annäherung an die zugrunde liegenden Daten, und im Wesentlichen schätzen Sie den Mittelwert und die Varianz auf Umwegen, indem Sie das Histogramm anstelle der Daten anpassen.

In jedem Fall, wenn Sie mit Ihrem Gedankengang weiter oben fortsetzen möchten, müssen Sie mehr Punkte zu den Winkeln hinzufügen. Der beste Weg, dies zu tun wäre, etwas zu tun, wie

angles2 = np.arange(-8,8,.1); 
plt.plot(angles2,gaus(angles2,*popt),'r',label='Fit') 

Es könnte sein, dass Ihre Passform sieht einfach schlecht aus, weil Sie sehr wenige Datenpunkte haben. Mit diesem Ansatz würden Sie sehen, wie die kontinuierliche Verteilung aussehen sollte.

Answer 2

Ich denke, es gibt zwei verschiedene Dinge hier:

scheinen eine Gaußsche Verteilung

→ zu folgen Wenn Sie denken, dass die Daten normal verteilt sind, Sie sind in den Bereichen Statistik und Wahrscheinlichkeitsverteilungen, und möchten vielleicht eine test machen, um zu sehen, ob sie mit einer bestimmten Verteilung (normal oder anders) übereinstimmen.

---

Und arbeiten mit Ihrem Grundstück:

erhalten ein "besseres" Gaußsche Grundstück

in Ihrem Code können Sie die erste Schätzung in curve_fit auslassen und Zeichnen Sie die angepasste Kurve gegen eine kontinuierliche unabhängige Variable:

import numpy as np 
import matplotlib.pyplot as plt 
from scipy import asarray as ar, exp, sqrt 
from scipy.optimize import curve_fit 


angles = [-8, -6, -4, -2, 0, 2, 4, 6, 8] 
data = [99, 610, 1271, 1804, 1823, 1346, 635, 125, 24] 
angles = ar(angles) 
data = ar(data) 

n = len(data) ## <--- 
mean = sum(data*angles)/n 
sigma = sqrt(sum(data*(angles-mean)**2)/n) 

def gaus(x,a,mu,sigma): 
    return a*exp(-(x-mu)**2/(2*sigma**2)) 

popt,pcov = curve_fit(gaus,angles,data)#,p0=[0.18,mean,sigma]) ## <--- leave out the first estimation of the parameters 
xx = np.linspace(-10, 10, 100) ## <--- calculate against a continuous variable 

fig = plt.figure() 
plt.plot(angles, data, "ob", label = "Measured") 
plt.plot(xx,gaus(xx,*popt),'r',label='Fit') ## <--- plot against the contious variable 
plt.xlim(-10, 10) 
plt.ylim(0, 2000) 
plt.xticks(angles) 
plt.title("$^{137}$Cs Zero Point") 
plt.xlabel("Angle [$^\circ$]") 
plt.ylabel("662 keV-Photon Count") 
plt.grid() 
plt.legend() 
plt.savefig('normal.png') 
plt.show() 

---

In diesem Beispiel:

print(popt) 

[ 1.93154077e+03 -9.21486804e-01 3.26251063e+00] 

anzumerken, dass die erste Schätzung des Parameters Grßenordnungen vom Ergebnis: 0.18 vs. 1931,15.