【问题标题】：Exponential decay curve fitting in numpy and scipynumpy和scipy中的指数衰减曲线拟合
【发布时间】：2012-12-30 06:17:20
【问题描述】：

我在将曲线拟合到某些数据时遇到了一些麻烦，但无法弄清楚我哪里出错了。

过去我使用 numpy.linalg.lstsq 来处理指数函数，使用 scipy.optimize.curve_fit 来处理 sigmoid 函数。这次我希望创建一个脚本，让我指定各种函数、确定参数并测试它们对数据的拟合度。在执行此操作时，我注意到 Scipy leastsq 和 Numpy lstsq 似乎为相同的数据集和相同的功能提供了不同的答案。该函数只是y = e^(l*x)，并被限制为y=1x=0。

Excel 趋势线与 Numpy lstsq 的结果一致，但由于 Scipy leastsq 能够接受任何功能，因此最好找出问题所在。

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

## Sampled data
x = np.array([0, 14, 37, 975, 2013, 2095, 2147])
y = np.array([1.0, 0.764317544, 0.647136491, 0.070803763, 0.003630962,     0.001485394,     0.000495131])

# function
fp = lambda p, x: np.exp(p*x)

# error function
e = lambda p, x, y: (fp(p, x) - y)

# using scipy least squares
l1, s =  optimize.leastsq(e, -0.004, args=(x,y))
print l1
# [-0.0132281]


# using numpy least squares
l2 = np.linalg.lstsq(np.vstack([x, np.zeros(len(x))]).T,np.log(y))[0][0]
print l2
# -0.00313461628963 (same answer as Excel trend line)

# smooth x for plotting
x_ = np.arange(0, x[-1], 0.2)

plt.figure()
plt.plot(x, y, 'rx', x_, fp(l1, x_), 'b-', x_, fp(l2, x_), 'g-')
plt.show()

编辑 - 附加信息

上面的 MWE 包括数据集的一个小样本。在拟合实际数据时，scipy.optimize.curve_fit 曲线的 R^2 为 0.82，而 numpy.linalg.lstsq 曲线与计算所得的曲线相同通过 Excel，R^2 为 0.41。

【问题讨论】：

标签： numpy scipy curve-fitting least-squares exponential

【解决方案1】：

您正在最小化不同的误差函数。

当你使用numpy.linalg.lstsq时，被最小化的误差函数是

np.sum((np.log(y) - p * x)**2)

scipy.optimize.leastsq 最小化函数

np.sum((y - np.exp(p * x))**2)

第一种情况需要因变量和自变量之间存在线性依赖关系，但解是通过分析已知的，而第二种情况可以处理任何依赖关系，但依赖于迭代方法。

另外说明，~~我现在无法测试它，但是~~当使用numpy.linalg.lstsq 时，我不需要vstack 一行零，以下也可以：

l2 = np.linalg.lstsq(x[:, None], np.log(y))[0][0]

【讨论】：

谢谢@Jaime - 很好的答案！不幸的是，我的数学知识不是那么好；是一个写还是错[也见上面的编辑]，还是他们只是根本不同......？对其他函数有什么影响，例如，如果我想测试 Sigmoid 或 Gompertz 曲线对相同数据的拟合？
@StacyR 我没有正确回答您的问题的知识，但我很确定像您对np.linalg.lstsq 所做的那样拟合指数只是一个快速的'n'dirty技巧，它不会不能正确计算错误。这里有一些讨论（我很难理解）：mathworld.wolfram.com/LeastSquaresFittingExponential.html 如果你不想深入研究这些东西，我会用 scipy 的方法来处理所有事情：它应该更适合，你的结果将是对所有功能保持一致。
再次感谢！我对此进行了更多研究，正如您所提到的，我发现np.linalg.lstsq 方法在低 x 值下过度重视 y 误差。您共享的链接以及我发现的其他一些资源使我能够推导出另一种分析方法（使其变得棘手的是约束 --- 所有书籍都描述了 y=ae^bx 而不是 y=e^b*x)，但是，这也会产生比迭代 scipy.optimize.leastsq 更差的拟合曲线。

【解决方案2】：

稍微解释一下 Jaime 的观点，数据的任何非线性变换都会导致不同的误差函数，从而导致不同的解决方案。这些将导致拟合参数的不同置信区间。因此，您可以使用三个可能的标准来做出决定：您想要最小化哪个错误，您希望对哪些参数更有信心，最后，如果您使用拟合来预测某个值，哪种方法在有趣的情况下产生的错误更少预测值。稍微分析一下并在 Excel 中进行分析表明，数据中不同种类的噪声（例如，如果噪声函数缩放幅度、影响时间常数或相加）会导致不同的解决方案选择。

我还要补充一点，虽然这个技巧“适用”于指数衰减到 0，但它不能用于更一般（和常见）的阻尼指数（上升或下降）到无法假设的值的情况为 0。

【讨论】：