在Python II中擬合高斯曲線

Question

我有兩個列表。

import numpy 
x = numpy.array([7250, ... list of 600 ints ... ,7849]) 
y = numpy.array([2.4*10**-16, ... list of 600 floats ... , 4.3*10**-16]) 

它們形成U形曲線。 現在我想適應高斯曲線。

from scipy.optimize import curve_fit 
n = len(x) 
mean = sum(y)/n 
sigma = sum(y - mean)**2/n 

def gaus(x,a,x0,sigma,c): 
    return a*numpy.exp(-(x-x0)**2/(2*sigma**2))+c 

popt, pcov = curve_fit(gaus,x,y,p0=[-1,mean,sigma,-5]) 

pylab.plot(x,y,'r-') 
pylab.plot(x,gaus(x,*popt),'k-') 
pylab.show() 

我剛結束了嘈雜的原始U形曲線和一條直線水平線穿過曲線。

我不確定-1和-5在上面的代碼中代表什麼，但我確定我需要調整它們或其他東西來獲得高斯曲線。我一直在玩弄可能的價值，但無濟於事。

任何想法？

Answer 1

p0在您調用curve_fit的過程中給出了除x之外的其他參數的初始猜測。在上面的代碼中，您要說curve_fit函數使用-1作爲a的初始猜測，-5作爲c的初始猜測，mean作爲x0的初始猜測，而sigma作爲sigma的猜測。然後curve_fit函數將調整這些參數以嘗試更好地擬合。問題在於你的函數參數的初始猜測在（x，y）s的順序上是非常糟糕的。

想想你的高斯的不同參數的數量級。 a應該在y值的大小（10 ** - 16）附近，就像高斯的峯值一樣，指數部分永遠不會大於1. x0會給出您x值的位置，在這個位置x值的指數部分高斯將是1，所以x0應該在7500左右，可能在數據的中心。西格瑪指出了你的高斯的寬度或擴散，所以也許是100年代的一些猜測。最後，c只是將整個高斯向上和向下移動的偏移量。

我推薦做的事情是在擬合曲線之前，爲a，x0，sigma和c選擇一些看似合理的值，並用高斯函數繪製數據，然後用a，x0，sigma和c，直到你得到的東西看起來至少有一些你希望高斯擬合的方式，然後用這些作爲curve_fit p0值的起點。我給的價值應該讓你開始，但可能不會完全按照你的意願去做。例如，如果要翻轉高斯以獲得「U」形狀，則可能需要爲負值。

另外打印出curve_fit認爲對您的a，x0，sigma和c有好處的值可能會幫助您瞭解它在做什麼以及該功能是否在正確的軌道上以最小化擬合殘差。

我有類似的問題做曲線擬合gnuplot，如果初始值太離你想要適應它完全錯誤的方向與參數來最大限度地減少殘差，你可以做得更好眼。把這些功能看作是通過眼睛估計這些參數來微調你的方法。

希望幫助

Answer 2

首先，你的變量sigma實際上是變化的，即西格瑪平方--- http://en.wikipedia.org/wiki/Variance#Definition。 通過給curve_fit一個次優的開始估計值，這會使其混淆。

然後，您的配件ansatz gaus包含振幅a和偏移量，這是您實際需要的嗎？起始值是a=-1（取消鐘形）並且偏移c=-5。他們來自哪裏？

這裏就是我想要做的：

• 解決您的擬合模型。你想只是一個高斯，它是否需要規範化。如果是這樣，則振幅a由sigma等修復。
• 查看實際數據。什麼是尾巴（偏移），什麼是符號（振幅符號）。

如果你真的想只是沒有任何花俏高斯，你可能實際上並不需要curve_fit：高斯由兩個第一時刻，mean和sigma完全定義。像你一樣計算它們，將它們繪製在數據上並查看你是否全部設置好。

Answer 3

感謝大家的解答。

x = numpy.array([7250, ... list of 600 ints ... ,7849]) 
y = numpy.array([2.4*10**-16, ... list of 600 floats ... , 4.3*10**-16]) 

n = len(x) 
mean = sum(x)/n 
sigma = math.sqrt(sum((x-mean)**2)/n) 

def gaus(x,a,x0,sigma,c): 
    return a*numpy.exp(-(x-x0)**2/(2*sigma**2))+c 

popy, pcov = curve_fit(gaus,x,y,p0=[-max(y),mean,sigma,min(x)+((max(x)-min(x)))/2]) 

pylab.plot(x,gaus(x,*popt)) 

Answer 4

我不認爲你是正確估計你的平均和西格瑪的初步猜測。

看看在SciPy的食譜here

我想應該是這樣的。

x = numpy.array([7250, ... list of 600 ints ... ,7849]) 
y = numpy.array([2.4*10**-16, ... list of 600 floats ... , 4.3*10**-16]) 

n = len(x) 
mean = sum(x*y)/sum(y) 
sigma = sqrt(abs(sum((x-mean)**2*y)/sum(y))) 

def gaus(x,a,x0,sigma,c): 
    return a*numpy.exp(-(x-x0)**2/(2*sigma**2))+c 

popy, pcov = curve_fit(gaus,x,y,p0=[-max(y),mean,sigma,min(x)+((max(x)-min(x)))/2]) 

pylab.plot(x,gaus(x,*popt)) 

如果任何人有一個簡單的解釋，爲什麼這些是正確的時刻的鏈接，我將不勝感激。我相信SciPy Cookbook是對的。