scipy.optimize.curvefit 使用边界时失败

Question

我正在尝试使用scipy.optimize.curvefit 将一组数据与一个函数（参见下面的示例）相匹配， 但是当我使用边界（documentation）时，拟合失败了，我只是得到 初始猜测参数作为输出。 一旦我将-np.inf 广告np.inf 替换为第二个参数的界限 （函数中的dt），适合的作品。 我做错了什么？

import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as opt

#Generate data
crc=np.array([-1.4e-14, 7.3e-14, 1.9e-13, 3.9e-13, 6.e-13, 8.0e-13, 9.2e-13, 9.9e-13, 
              1.e-12, 1.e-12, 1.e-12, 1.0e-12, 1.1e-12, 1.1e-12, 1.1e-12, 1.0e-12, 1.1e-12])
time=np.array([0., 368., 648., 960., 1520.,1864., 2248., 2655., 3031., 
                3384., 3688., 4048., 4680., 5343., 6055.,  6928.,  8120.])

#Define the function for the fit
def testcurve(x, Dp, dt):
    k = -Dp*(x+dt)*2e11
    curve = 1e-12 * (1+2*(-np.exp(k) + np.exp(4*k) - np.exp(9*k) + np.exp(16*k)))
    curve[0]= 0
    return curve

#Set fit bounds 
dtmax=time[2]
param_bounds = ((-np.inf, -dtmax),(np.inf, dtmax))

#Perform fit
(par, par_cov) = opt.curve_fit(testcurve, time, crc, p0 = (5e-15, 0), bounds = param_bounds)

#Print and plot output
print(par)       
plt.plot(time, crc, 'o')
plt.plot(time, testcurve(time, par[0], par[1]), 'r-')
plt.show()

Answer 1

我今天在不同的拟合问题中遇到了相同的行为。在网上搜索了一下，发现这个链接很有帮助：Why does scipy.optimize.curve_fit not fit to the data?

简短的回答是：在数值拟合中使用极小（或大）的数字并不可靠，缩放它们会导致更好的拟合。

---

在您的情况下，crc 和 Dp 都是可以放大的极小数字。您可以使用比例因子，并且在一定范围内，拟合看起来非常健壮。完整示例：

import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as opt

#Generate data
crc=np.array([-1.4e-14, 7.3e-14, 1.9e-13, 3.9e-13, 6.e-13, 8.0e-13, 9.2e-13, 9.9e-13, 
              1.e-12, 1.e-12, 1.e-12, 1.0e-12, 1.1e-12, 1.1e-12, 1.1e-12, 1.0e-12, 1.1e-12])
time=np.array([0., 368., 648., 960., 1520.,1864., 2248., 2655., 3031., 
                3384., 3688., 4048., 4680., 5343., 6055.,  6928.,  8120.])

# add scale factors to the data as well as the fitting parameter
scale_factor_1 = 1e12 # 1./np.mean(crc) also works if you don't want to set the scale factor manually
scale_factor_2 = 1./2e11

#Define the function for the fit
def testcurve(x, Dp, dt):
    k = -Dp*(x+dt)*2e11 * scale_factor_2
    curve = 1e-12 * (1+2*(-np.exp(k) + np.exp(4*k) - np.exp(9*k) + np.exp(16*k))) * scale_factor_1
    curve[0]= 0
    return curve

#Set fit bounds 
dtmax=time[2]
param_bounds = ((-np.inf, -dtmax),(np.inf, dtmax))

#Perform fit
(par, par_cov) = opt.curve_fit(testcurve, time, crc*scale_factor_1, p0 = (5e-15/scale_factor_2, 0), bounds = param_bounds)

#Print and plot output
print(par[0]*scale_factor_2, par[1])

plt.plot(time, crc*scale_factor_1, 'o')
plt.plot(time, testcurve(time, par[0], par[1]), 'r-')
plt.show()

拟合结果：[6.273102923176595e-15, -21.12202697564494]，它给出了一个合理的拟合，并且非常接近没有任何界限的结果：[6.27312512e-15, -2.11307470e+01]