3

我正在寻找在拟合参数中输出不确定性的最简单方法。使用 spo.curve_fit,我们只需在拟合时得到协方差矩阵,我们可以取对角线和平方根来找出不确定性。使用 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 打印参数的不确定性,就像上面最后一行代码一样,但这只是返回“无”。我可以得到一个协方差矩阵,但我不知道它是以什么顺序显示的。我的最终目标只是获得值和相关的不确定性,然后我可以将它们放入一个数组中并在我的代码中进一步使用。

4

2 回答 2

2

不确定性隐藏在 .stderr 属性下的模型中。例如 :

gmodel = Model(function_tofit)
fit_params = Parameters()
fit_params.add('parameter1',value=guess1)
fit_params.add('parameter2',value=guess2)
...
result = gmodel.fit(y_to_fit, fit_params, x = x_to_fit)

#get the output value of the fit
mean_fit_value = result.params['parameter1'].value
std_fit_value = result.params['parameter1'].stderr

在您的情况下,您打印了 fit_params['your_parameter'].stderr 这显然是 None 因为您没有在参数上指定任何先验。请记住 fit_params 是您的输入,而 params 是您想要的输出。

于 2019-09-02T21:54:55.140 回答
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...

干杯!

于 2015-05-19T05:30:21.960 回答