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

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 に答える