5

这些天我花了一些时间在一个问题上。我有一组数据：

y = f(t)，其中 y 是非常小的浓度 (10^-7)，t 排在第二位。t 从 0 到 12000 左右变化。

测量遵循既定模型：

y = Vs * t - ((Vs - Vi) * (1 - np.exp(-k * t)) / k)

我需要找到 Vs、Vi 和 k。所以我使用了curve_fit，它返回了最佳拟合参数，并绘制了曲线。

然后我使用了一个类似的模型：

y = (Vs * t/3600 - ((Vs - Vi) * (1 - np.exp(-k * t/3600)) / k)) * 10**7

通过这样做，t 是小时数，y 是介于 0 和大约 10 之间的数字。返回的参数当然是不同的。但是当我绘制每条曲线时，我得到的是：

http://i.imgur.com/XLa4LtL.png

绿色拟合是第一个模型，蓝色拟合是“标准化”模型。红点是实验值。

拟合曲线不同。我认为这不是预期的，我不明白为什么。如果数字“合理”，计算是否更准确？

4

3 回答 3

5

的文档字符串optimize.curve_fit说，

p0 : None, scalar, or M-length sequence
    Initial guess for the parameters.  If None, then the initial
    values will all be 1 (if the number of parameters for the function
    can be determined using introspection, otherwise a ValueError
    is raised).

因此，首先，参数的初始猜测默认为 1。

此外，曲线拟合算法必须针对各种参数值对函数进行采样。最初选择“各种值”时的初始步长约为 1。如果您的数据随着参数值的变化而平滑地变化（大约为 1），则该算法将运行得更好。

如果函数随着 1 级的参数变化而剧烈变化，那么算法可能会倾向于错过最佳参数值。

请注意，即使该算法在调整参数值时使用了自适应步长，如果初始调整离目标太远而产生较大的残差，并且如果在其他方向上的调整恰好产生较小的残差，那么该算法可能会在错误的方向上徘徊并错过局部最小值。它可能会找到其他一些（不需要的）局部最小值，或者根本无法收敛。因此，使用具有自适应步长的算法不一定会节省您的时间。

这个故事的寓意是，扩展数据可以提高算法找到所需最小值的机会。

当应用于数量级为 1 的数据时，数值算法通常都倾向于更好地工作。这种偏差以多种方式进入算法。例如，optimize.curve_fit依赖optimize.leastsq，其调用签名为optimize.leastsq：

def leastsq(func, x0, args=(), Dfun=None, full_output=0,
            col_deriv=0, ftol=1.49012e-8, xtol=1.49012e-8,
            gtol=0.0, maxfev=0, epsfcn=None, factor=100, diag=None):

因此，默认情况下，公差ftol和xtol的数量级为 1e-8。如果找到最佳参数值需要小得多的容差，那么这些硬编码的默认数字将导致optimize.curve_fit错过优化参数值。

为了使这一点更具体，假设您试图最小化f(x) = 1e-100*x**2. 1e-100 的因子将 - 值压缩得y如此之多，以至于大范围的 -x值（上面提到的参数值）将适合 1e-8 的容差。因此，如果缩放不理想，leastsq将无法很好地找到最小值。

使用大约 1 的浮点数的另一个原因是，间隔中的 (IEEE754) 浮点数[-1,1]比远离 1 的浮点数多得多。例如，

import struct
def floats_between(x, y):
    """
    http://stackoverflow.com/a/3587987/190597 (jsbueno)
    """
    a = struct.pack("<dd", x, y)
    b = struct.unpack("<qq", a)
    return b[1] - b[0]

In [26]: floats_between(0,1) / float(floats_between(1e6,1e7))
Out[26]: 311.4397707054894

这表明表示 0 和 1 之间数字的浮点数是区间 [1e6, 1e7] 中的浮点数的 300 多倍。因此，在其他条件相同的情况下，如果使用小数字而不是非常大的数字，您通常会得到更准确的答案。

于 2013-09-06T17:49:19.427 回答
2

我想这与您传递给曲线拟合的初始参数估计有关。如果你没有通过任何我相信他们都默认为 1。标准化你的数据使这些初始估计更接近真相。如果您不想使用标准化数据，只需自己传递初始估计值并给它们合理的值。

于 2013-09-06T17:45:27.427 回答
2

其他人已经提到，您可能需要对自己的健康状况进行良好的初步猜测。在这种情况下，我通常会尝试找到一些快速而肮脏的技巧来至少对参数进行大致估计。在你的情况下，对于 large t，指数很快衰减到零，所以对于 large t，你有

y == Vs * t - (Vs - Vi)  / k

进行一阶线性拟合，例如

[slope1, offset1] = polyfit(t[t > 2000], y[t > 2000], 1)

你会得到slope1 == Vs和offset1 == (Vi - Vs) / k。

从你拥有的所有点中减去这条直线，你得到指数

residual == y - slope1 * t - offset1 == (Vs - Vi) * exp(-t * k)

取双方的日志，你得到

log(residual) == log(Vs - Vi) - t * k

所以做第二次合身

[slope2, offset2] = polyfit(t, log(y - slope1 * t - offset1), 1)

会给你slope2 == -kand offset2 == log(Vs - Vi)，这应该是可以解决的，Vi因为你已经知道了Vs。您可能必须将第二次拟合限制为 的较小值t，否则您可能会获取负数的对数。收集您通过这些拟合获得的所有参数，并将它们用作您的curve_fit.

最后，您可能想考虑进行某种加权拟合。关于曲线指数部分的信息仅包含在前几个点中，所以也许你应该给它们更高的权重。以统计上正确的方式做到这一点并非易事。

于 2013-09-06T18:11:16.693 回答