python - 将数据拟合到后处理非线性模型时，curve_fit 的算法是什么？

Question

我对编程工具和Python都很陌生，对数值计算知之甚少。我试图理解 Curve_fit 的源代码，但这对我来说太难了。

如果我们直接对非线性模型进行曲线拟合，通过计算 Hessian 矩阵，我想我们能够找到所有参数的临界点，无论局部或全局、最小值或最大值。

但是，如果包含模型的函数无法微分，那么 curve_fit 实际上如何计算答案？

我想象该算法正在寻找所有数据点的每个参数的每个区间中 Chi^2 的局部最小值。如果这样做，默认间隔有多大，每次参数试验之间的距离是多少，最大迭代是多少？如果像我说的那样，那么这个迭代在多个参数之间是如何工作的呢？他们每次都是作为一组单独尝试吗？

我想了解该算法，以便我可以编写一个好的函数来应用 curve_fit。我试图适应的函数中存在许多脏修复，不同的小更改中出现不同的错误，因此我无法将代码放在这里，因为我不知道问题所在。

另外，关于 sigma 输入，如果我没有 y 错误的线索，如果默认 sigma 与数据值相比太大，结果会怎样？或者，如果模型的灵 active 对后处理影响不大，或者 sigma 远小于数据到函数的距离，会发生什么？

Answer 1

要开始了解 scipy curve_fit 如何尝试解决问题，请先阅读非线性优化方法，例如 Levenberg-Marquardt (LM) 方法(请参阅 https://en.wikipedia.org/wiki/Levenberg%E2%80%93Marquardt_algorithm )。这是未指定参数值界限时 curve_fit 使用的方法，它将给出您寻求的大部分答案。

此方法(以及 scipy.optimize 中的其他求解器)从初始值开始，并尝试细化这些值，以找到最小结果残差数组(在最小二乘意义上:需要明确的是，该方法需要一个残差数组并自行进行平方和求和)。为了确定参数值移动的方向，它使用一阶导数(变量的残差)或雅可比行列式。对此没有分析形式是很常见的，在这种情况下，通过在参数值中采取非常小的步长(通常围绕机器精度)来计算数值导数。

利用这些导数和最后的步长大小，LM 方法决定需要多大的步长才能找到最小残差。将残差视为每个参数的二次函数是一种常见的(而且通常是好的)近似。如果该值远离底部，则采取导数呈线性的步骤(即，仅遵循局部斜率)就非常好。但在接近最终解决方案时，使用二阶导数(又名 Hessian、曲率或某种意义上的“历史”)会很有帮助。 LM 结合这两种方法来尝试快速找到解决方案。在每个步骤中，它首先计算导数(以及第一步之后的二阶导数)，然后计算所有参数要采取的步骤。它重复此过程，直到解(卡方)与容差相比没有改善。在scipy的leastsq(通常由curve_fit使用)中，可以指定用于计算数值导数的步长以及停止拟合的容差。

通常，当使用 curve_fit 或 leastsq 或 least_squares 时，您无需担心任何有关解决方案如何实现的细节成立。好吧，除了如果您确实有分析导数，它可以提高过程的速度和准确性。您的“模型函数”应该只接受传入的参数并计算数据的预测模型。更改传入的值会使算法感到困惑 - 它会使用它知道的值运行您的函数，以尝试找到解决方案。

其他问题:
如果您不知道数据中的标准错误 (y)，请不要太担心。报告卡方的绝对值可能不会接近 ndata-nvarys(正如预期的良好拟合)，但相对于其他拟合的值应该很好。

LM 方法的特征之一(至少在 MINPACK -> leastsq -> curve_fit 实现中)是它报告最终的协方差矩阵，该矩阵可以用于确定拟合参数的不确定性和相关性。通常，报告的值应将卡方增加 1(即 1-sigma 标准误差)。由于数据缺乏良好的不确定性非常常见，因此 scipy curve_fit 会缩放此协方差以给出值，就好像卡方等于 ndata-nvarys 一样。也就是说，如果您认为拟合良好，则参数中的标准误差将相当不错(并且可以使用 sqrt(chi-square/(ndata-nvays)) 来估计参数中的误差数据)。