我正在使用python lmfit包进行数据拟合。我正在观察几个峰值,并将它们与伪Voigt+常数模型拟合。然后,我提取拟 symfony 的中心位置,并编译大约一百个这样的峰位置。然后用正弦波模型(也用lmfit)绘制和拟合一百个位置。我的问题是,我无法找到显示从伪Voigt拟合提取的峰值位置的不确定性所需的信息。在许多类似的问题中,我遇到了同样包含在lmfit包中的eval_determination函数,但我无法使用它来提取中心位置的置信区间,只能提取峰值的振幅和σ值的不确定性。
理想的结果如下所示:The ideal case I'm trying to achieve
我现在的图表看起来是这样的:Current graph的
显然,我需要从散点图切换到误差条图,但我目前没有任何误差区间可供使用。
以下是我能构建的最简单可行的“迷你”案例,同时仍然包括我原来的适合模型。我的原始数据大约是100个数据点,而不是4个,因此在这个较小的情况下,正弦波不是“平滑”的,但这不应该影响引入误差条的方法。
import numpy as np
import matplotlib.pyplot as plt
from lmfit.models import PseudoVoigtModel, ConstantModel, SineModel
def Voigt_app(x,y,a,b):
psvo_mod = PseudoVoigtModel(prefix='psvo_')
pars = psvo_mod.guess(y,x=x)
con_mod = ConstantModel(prefix='con_')
pars = pars.update(con_mod.make_params(c = (y[0]+y[1])/2))
mod = psvo_mod + con_mod
out = mod.fit(y, pars, x=x)
details = out.fit_report(min_correl=0.25)
return out, details
def Sin_app(x,y):
sin_mod = SineModel(prefix='sin_')
pars = sin_mod.guess(y,x=x)
con_mod = ConstantModel(prefix='con_')
pars = pars.update(con_mod.make_params(c = np.mean(y)))
mod = sin_mod + con_mod
out = mod.fit(y,pars,x=x)
details = out.fit_report(min_correl=0.25)
return out, details
data = np.asarray(
[[14407.0 , 26989.0 , 31026.0 , 31194.0 , 31302.0 , 91996.0 , 112250.0 , 112204.0 , 113431.0 , 75097.0 , 55942.0 , 55826.0 , 56181.0 , 28266.0 , 14687.0],
[12842.0 , 13275.0 , 13581.0 , 13943.0 , 14914.0 , 20463.0 , 21132.0 , 21457.0 , 23308.0 , 32017.0 , 30808.0 , 30927.0 , 30337.0 , 21025.0 , 15216.0],
[15770.0 , 17677.0 , 19008.0 , 20911.0 , 25958.0 , 27984.0 , 28164.0 , 31161.0 , 33517.0 , 29430.0 , 28673.0 , 32725.0 , 28495.0 , 24527.0 , 25173.0],
[16299.0 , 20067.0 , 25102.0 , 34968.0 , 37757.0 , 44871.0 , 51347.0 , 60154.0 , 54985.0 , 53383.0 , 45776.0 , 40836.0 , 30816.0 , 27922.0 , 26066.0]])
a,b = 1932, 1947
centers = []
colors = ['red','black','green','blue']
fig, ax = plt.subplots()
for i in range(4):
ax.scatter(np.arange(15),data[i,:],label=f'data {i}',color=colors[i])
out, details = Voigt_app(x,data[i,:],a,b)
ax.plot(out.best_fit,label=f'fit {i}',color=colors[i])
centers.append(out.values['psvo_center'])
ax.set_title('Data points and pseudo-Voigt fits')
ax.legend()
fig, ax = plt.subplots()
ax.scatter(np.arange(4),centers,color='red',label='data')
out, details = Sin_app(np.arange(4),np.asarray(centers))
ax.plot(np.arange(4),out.best_fit,color='black',label='fit')
ax.set_title('Peak positions and sine-wave fit')
ax.legend()
字符串
由此,我想用误差条形图切换第二个图中的散点图,其中包含伪Voigt峰位置的置信区间(“out.values['psvo_center']”)。
2条答案
按热度按时间xlpyo6sf1#
如果我对你的问题理解正确的话,你不仅要提取每次拟合的'psvo_center'参数的最佳拟合值,还要提取'psvo_center'的标准误差,然后使用这些数据来拟合正弦振荡。
如果这是正确的,我认为您需要做的是保存
stderr
属性,方法是将代码改为字符串
然后,您可以使用它来制作和误差条图
型
您可能还希望使用这些
stderrs
来对正弦函数的拟合进行加权,也许可以将Sin_App
修改为型
然后用
型
对于你给予的4个示例数据集,标准误差的变化不是很大,所以它不会产生很大的差异。但随着您添加更多数据集,它可能会变得更加重要。您可能希望在使用此“in production”之前添加一个检查,以确保
stderr
有效(例如,不是None
,也不是0
)。对于“另一件事”,您可能还希望计算正弦函数中的预测不确定性,并绘制该范围。您的示例中有4个正弦拟合变量,但只有4个值,因此在这种情况下,拟合正弦函数的不确定性非常小,几乎完美。因此,为了说明该方法,我在下面的拟合中修复了
con_c
参数。如果拟合超过4个峰形的结果,则不需要。但是,为了说明,将Sin_app
更改为型
然后把结果和
型
这导致最终图:
x1c 0d1x的数据
示出了具有误差条的峰中心的“数据”、拟合值和拟合的1-σ范围。
cngwdvgl2#
以下是如何修改Voigt_app函数来实现这一点:
字符串
然后,您可以在绘制误差线时使用这些不确定性:
型