Scipy 曲线拟合结果异常

Question

当我尝试拟合我的数据时，结果有点奇怪，我不明白为什么。得到的拟合结果很平坦，而且第一个输入 e=0. 似乎在某个地方引发了除零错误。唯一有效的情况是当我把 e[0] 修改为 1.0e-9。

结果如下所示：

从这个 例子 来看，我的例子似乎和我读到的内容差不多，但我还是卡住了，所以你能帮我看看我的情况出错在哪里吗？

import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

src_s = np.array((45.59,50.66664,59.74871,65.71018,72.76012,79.06256,84.13755,90.39944,
                  96.33653,101.65667,106.27968,110.76301,114.41808,117.21922,120.51836))
src_e = np.array((0.0,0.00126,0.00503,0.00804,0.01228,0.01685,0.02127,0.02846,0.03666,
                  0.04581,0.05620,0.06882,0.08005,0.09031,0.10327))
# plot source data
plt.plot(src_e, src_s, 'o')
# fitting function
def sigma(e, k ,n): return k*(e**n)
# find parameters curve fitting
param, var = curve_fit(sigma, src_e, src_s)
new_e = np.linspace(src_e.min(), src_e.max(), 50)
plt.plot(new_e, sigma(new_e, *param))

# modify first input
src_e[0]=1.0e-9
# relaunch parameters curve fitting
param, var = curve_fit(sigma, src_e, src_s)
new_e = np.linspace(src_e.min(), src_e.max(), 50)
plt.plot(new_e, sigma(new_e, *param))

plt.show()

提前感谢你的帮助。

Answer 1

第一个点不能在曲线上，所以你需要修改曲线的公式：

import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

src_s = np.array((45.59,50.66664,59.74871,65.71018,72.76012,79.06256,84.13755,90.39944,
                  96.33653,101.65667,106.27968,110.76301,114.41808,117.21922,120.51836))
src_e = np.array((0.0,0.00126,0.00503,0.00804,0.01228,0.01685,0.02127,0.02846,0.03666,
                  0.04581,0.05620,0.06882,0.08005,0.09031,0.10327))
# plot source data
plt.plot(src_e, src_s, 'o')

def sigma(e, k ,n, offset): return k*((e+offset)**n)
# find parameters curve fitting
param, var = curve_fit(sigma, src_e, src_s)
new_e = np.linspace(src_e.min(), src_e.max(), 50)
plt.plot(new_e, sigma(new_e, *param))

这是输出结果：

Answer 2

问题的根本原因是参数的初始猜测不太好（实际上没有为 curve_fit 提供起始参数）。

我们的目标函数可以很容易地进行线性化处理。我们可以先把它线性化，然后进行线性回归，这样就能得到一组比较好的初始参数，传给 curve_fit（通过 p0= 传递）。这样得到的拟合效果更好（残差更小），而且不需要把第一个值替换成 1e-9：

In [38]:

src_e[0]=1.0e-9
# relaunch parameters curve fitting
param, var = curve_fit(sigma, src_e, src_s)
new_e = np.linspace(src_e.min(), src_e.max(), 50)
src_e[0]=0
plt.plot(new_e, sigma(new_e, *param))
plt.plot(src_e, src_s, 'ro')
plt.savefig('1.png')
print 'Residue is:', ((sigma(src_e, *param)-src_s)**2).sum()
Residue is: 2168.65307587

In [39]:

import scipy.stats as ss
src_e[0]=0
V=ss.linregress(np.log(src_e)[1:], np.log(src_s)[1:]) #avoid log(0)
param, var = curve_fit(sigma, src_e, src_s, p0=(np.exp(V[1]), V[0]))
new_e = np.linspace(src_e.min(), src_e.max(), 50)
plt.plot(new_e, sigma(new_e, *param))
plt.plot(src_e, src_s, 'ro')
plt.savefig('1.png')
print 'Residue is:', ((sigma(src_e, *param)-src_s)**2).sum()
Residue is: 2128.85364181