python

2024-03-16 15:53| 来源: 网络整理| 查看: 265

我正在尝试在 python 中拟合一些具有 y 不确定性的数据点。数据在 python 中标记为 x、y 和 yerr。我需要以 loglog 比例对该数据进行线性拟合。作为拟合结果是否正确的引用，我将 python 结果与 Scidavis 的结果进行了比较

我试过 curve_fit

def func(x, a, b): return np.exp(a* np.log(x)+np.log(b)) popt, pcov = curve_fit(func, x, y,sigma=yerr)

以及 kmpfit 与

def funcL(p, x): a,b = p return ( np.exp(a*np.log(x)+np.log(b)) ) def residualsL(p, data): a,b=p x, y, errorfit = data return (y-funcL(p,x)) / errorfit a0=1 b0=0.1 p0 = [a0,b0] fitterL = kmpfit.Fitter(residuals=residualsL, data=(x,y,yerr)) fitterL.parinfo = [{}, {}] fitterL.fit(params0=p0)

当我试图用其中一个没有不确定性的数据(即设置 yerr=1)拟合数据时，一切正常，结果与 scidavis 的结果相同。但是如果我将 yerr 设置为数据文件的不确定性，我会得到一些令人不安的结果。在 python 中，我得到即 a=0.86，在 scidavis 中，a=0.14。我读到一些关于错误包含在权重中的信息。我是否必须更改任何内容才能正确计算拟合度？或者我做错了什么？

编辑:这里是数据文件 (x,y,yerr) 的示例

3.942387e-02 1.987800e+00 5.513165e-01 6.623142e-02 7.126161e+00 1.425232e+00 9.348280e-02 1.238530e+01 1.536208e+00 1.353088e-01 1.090471e+01 7.829126e-01 2.028446e-01 1.023087e+01 3.839575e-01 3.058446e-01 8.403626e+00 1.756866e-01 4.584524e-01 7.345275e+00 8.442288e-02 6.879677e-01 6.128521e+00 3.847194e-02 1.032592e+00 5.359025e+00 1.837428e-02 1.549152e+00 5.380514e+00 1.007010e-02 2.323985e+00 6.404229e+00 6.534108e-03 3.355974e+00 9.489101e+00 6.342546e-03 4.384128e+00 1.497998e+01 2.273233e-02

结果:

in python: without uncertainties: a=0.06216 +/- 0.00650 ; b=8.53594 +/- 1.13985 with uncertainties: a=0.86051 +/- 0.01640 ; b=3.38081 +/- 0.22667 in scidavis: without uncertainties: a = 0.06216 +/- 0.08060; b = 8.53594 +/- 1.06763 with uncertainties: a = 0.14154 +/- 0.005731; b = 7.38213 +/- 2.13653

最佳答案

我一定是误会了什么。您发布的数据看起来不像

f(x,a,b) = np.exp(a*np.log(x)+np.log(b))

红线是scipy.optimize.curve_fit的结果, 绿线是 scidavis 的结果。

我的猜测是这两种算法都没有收敛到一个很好的拟合，所以结果不匹配也就不足为奇了。

我无法解释 scidavis 如何找到它的参数，但根据我理解的定义，scipy正在寻找最小二乘残差低于 scidavis 的参数:

import numpy as np import matplotlib.pyplot as plt import scipy.optimize as optimize def func(x, a, b): return np.exp(a* np.log(x)+np.log(b)) def sum_square(residuals): return (residuals**2).sum() def residuals(p, x, y, sigma): return 1.0/sigma*(y - func(x, *p)) data = np.loadtxt('test.dat').reshape((-1,3)) x, y, yerr = np.rollaxis(data, axis = 1) sigma = yerr popt, pcov = optimize.curve_fit(func, x, y, sigma = sigma, maxfev = 10000) print('popt: {p}'.format(p = popt)) scidavis = (0.14154, 7.38213) print('scidavis: {p}'.format(p = scidavis)) print('''\ sum of squares for scipy: {sp} sum of squares for scidavis: {d} '''.format( sp = sum_square(residuals(popt, x = x, y = y, sigma = sigma)), d = sum_square(residuals(scidavis, x = x, y = y, sigma = sigma)) )) plt.plot(x, y, 'bo', x, func(x,*popt), 'r-', x, func(x, *scidavis), 'g-') plt.errorbar(x, y, yerr) plt.show()

产量

popt: [ 0.86051258 3.38081125] scidavis: (0.14154, 7.38213) sum of squares for scipy: 53249.9915654 sum of squares for scidavis: 239654.84276

关于python - 如何正确地包含与 python 拟合的不确定性，我们在Stack Overflow上找到一个类似的问题： https://stackoverflow.com/questions/11833538/

【本文地址】

python

python

今日新闻

推荐新闻