2015-05-07 80 views
1

我正在尋找輸出擬合參數中不確定度的最簡單方法。使用spo.curve_fit,當我們擬合時,我們得到協方差矩陣,我們可以用對角線和平方根來找出不確定性。 lmfit似乎並不那麼簡單。lmfit中擬合參數的不確定性

我的裝修看起來是這樣的:

import lmfit 


a_lm2 = lmfit.Parameter('a', value=a_est) 
b_lm2 = lmfit.Parameter('b', value=b_est) 
x0_core_lm2 = lmfit.Parameter('x0_core', value=gaus1['x0_core']) 
x0_1_lm2 = lmfit.Parameter('x0_1', value=gaus1['x0_1']) 
x0_2_lm2 = lmfit.Parameter('x0_2', value=gaus1['x0_2']) 
x0_3_lm2 = lmfit.Parameter('x0_3', value=gaus1['x0_3']) 
x0_4_lm2 = lmfit.Parameter('x0_4', value=gaus1['x0_4']) 
sig_core_lm2 = lmfit.Parameter('sig_core', value=gaus1['sig_core']) 
sig_1_lm2 = lmfit.Parameter('sig_1', value=gaus1['sig_1']) 
sig_2_lm2 = lmfit.Parameter('sig_2', value=gaus1['sig_2']) 
sig_3_lm2 = lmfit.Parameter('sig_3', value=gaus1['sig_3']) 
sig_4_lm2 = lmfit.Parameter('sig_4', value=gaus1['sig_4']) 
m_lm2 = lmfit.Parameter('m', value=m, vary=False) 
c_lm2 = lmfit.Parameter('c', value=c, vary=False) 


gausfit2 = mod.fit(y, x=x, a=a_lm2, b=b_lm2, x0_core=x0_core_lm2, x0_1=x0_1_lm2, x0_2=x0_2_lm2, 

x0_3=x0_3_lm2, x0_4=x0_4_lm2, sig_core=sig_core_lm2, sig_1=sig_1_lm2, sig_2=sig_2_lm2, 

sig_3=sig_3_lm2, sig_4=sig_4_lm2, m=m_lm2, c=c_lm2,weights=None, scale_covar=False) 


print 'a_lm2_unc =', a_lm2.stderr 

當我生成一個合適的報告,我得到的不確定性值,所以他們顯然被計算。我的問題是呼喚他們並使用他們。我試圖用上面代碼的最後一行打印stderr參數的不確定性,但這只是返回'None'。我可以得到一個協方差矩陣,但我不知道顯示的順序是什麼。我的最終目標是簡單地將值和相關的不確定性放入數組中,並在我的代碼中進一步使用。

回答

0

因爲我沒有你的數據,我無法測試它,但它看起來像你幾乎在那裏。您的「gausfit2」應該是ModelFit對象(http://cars9.uchicago.edu/software/python/lmfit/model.html#model.ModelFit)。因此,您需要做的生成報告如下:

print gausfit2.fit_report #will print you a fit report from your model 

#you should also be able to access the best fit parameters and other associated attributes. For example, you can use gausfit2.best_values to obtain a dictionary of the best fit values for your params, or gausfit2.covar to obtain the covariance matrix you are interested in. 

print gausfit2.covar 

#one suggestion to shorten your writing is to just create a parameters class with a shorter name and add your params to that. 

params = lmfit.Parameters() 
params.add('a', value=a_est) #and so on... 

乾杯!

+0

我忘了提及你檢查lmfit中的F-statistics選項(http://cars9.uchicago.edu/software/python/lmfit/confidence.html),因爲你的信心限制並不總是是對稱的(儘管在這個高斯的例子中它們可能是)。我會建議選擇一個更合適的F-stat,但是,由於內置會導致相對較小的置信限制。這裏是一個爲什麼這是重要的鏈接:http://stackoverflow.com/questions/30276163/how-to-understand-f-test-based-lmfit-confidence-intervals/30282340#30282340 – XtremeJake