可以用scipy拟合曲线最小化，但不能用scipy curve_fit

Question

我正在尝试使用y= 1-a(1-bx)**n将函数scipy curve_fit拟合到一些实验数据。该模型仅在y> 0时存在，因此我裁剪了计算值以强制执行此操作。代码如下所示

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

# Driver function for scipy.minimize

def driver_func(x, xobs, yobs):

    # Evaluate the fit function with the current parameter estimates

    ynew = myfunc(xobs, *x)
    yerr = np.sum((ynew - yobs) ** 2)

    return yerr

# Define function

def myfunc(x, a, b, n):

    y = 1.0 - a * np.power(1.0 - b * x, n) 
    y = np.clip(y, 0.00, None )

    return y

if __name__ == "__main__":

    # Initialise data

    yobs = np.array([0.005, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.004, 
                    0.048, 0.119, 0.199, 0.277, 0.346, 0.395, 0.444, 0.469, 
                    0.502, 0.527, 0.553, 0.582, 0.595, 0.603, 0.612, 0.599])
    xobs = np.array([0.013, 0.088, 0.159, 0.230, 0.292, 0.362, 0.419, 0.471,
                    0.528, 0.585, 0.639, 0.687, 0.726, 0.772, 0.814, 0.854,
                    0.889, 0.924, 0.958, 0.989, 1.015, 1.045, 1.076, 1.078])

    # Initial guess

    p0 = [2.0, 0.5, 2.0]

    # Check fit pre-regression

    yold = myfunc(xobs, *p0)
    plt.plot(xobs, yobs, 'ko', label='data', fillstyle='none')
    plt.plot(xobs, yold, 'g-', label='pre-fit: a=%4.2f, b=%4.2f, n=%4.2f' % tuple(p0))

    # Fit curve using SCIPY CURVE_FIT

    try:
        popt, pcov = scipy.optimize.curve_fit(myfunc, xobs, yobs, p0=p0)
    except:
        print("Could not fit data using SCIPY curve_fit")
    else:
        ynew = myfunc(xobs, *popt)
        plt.plot(xobs, ynew, 'r-', label='post-curve_fit: a=%4.2f, b=%4.2f, n=%4.2f' % tuple(popt))

    # Fit curve using SCIPY MINIMIZE

    res = scipy.optimize.minimize(driver_func, p0, args=(xobs, yobs), method='Nelder-Mead')
    ynw2 = myfunc(xobs, *res.x)
    plt.plot(xobs, ynw2, 'y-', label='post-minimize: a=%4.2f, b=%4.2f, n=%4.2f' % tuple(res.x))

    plt.legend()
    plt.show()

我还使用SCIPY MINIMIZE实现相同的目的。如下图所示，MINIMIZE可以工作，但是CURVE_FIT基本上用尽了精力，甚至放弃了，尽管开始的猜测与MINIMIZE解决方案相差不远（至少在视觉上）。感谢您对为何curve_fit似乎在这里不起作用的任何想法。

谢谢！

Answer 1

此问题是由函数定义中的裁剪引起的。两种最小化方法的工作原理根本不同，因此，对削波的反应也大不相同。这里minimize与Nelder-Mead一起使用，这是一种无梯度方法。因此，该算法不计算数值梯度，也不估计任何雅可比行列式。与之相反，最终由least-squares调用的curve_fit正是这样做的。但是，如果函数不连续，则近似梯度并由此得出任何雅可比行列式都存在一定问题。如前所述，这种不连续性是由np.clip引入的。删除后，可以很容易地看出P0的猜测不如包含裁剪的结果好。 curve_fit确实会随着maxfev=5000的增加而收敛，而minimize在将方法更改为method='CG'时会立即失败。要查看算法困难，可以尝试手动提供jac。