如何使用 scipy.log(a-x) 类型函数拟合 log(a-x) 类型的函数。曲线拟合？答案

【问题标题】：How to fit log(a-x) type functions with scipy. curve_fit?如何使用 scipy.log(a-x) 类型函数拟合 log(a-x) 类型的函数。曲线拟合？
【发布时间】：2018-04-26 04:45:22
【问题描述】：

我正在尝试拟合一个看起来像log(y)=a*log(b-x)+c 的函数，其中a、b 和c 是需要拟合的参数。相关的代码是

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

def logfunc(T, a, b, c):
    v=(a*np.log(b-T))+c
    return v

popt, pcov=curve_fit(logfunc, T, np.log(Energy), check_finite=False, bounds=([0.1, 1.8, 0.1], [1.0, 2.6, 1.0]))

plt.plot(T, logfunc(T, *popt))
plt.show

T 和 Energy 是生成的一些数据（我用它来绘制其他东西，所以数据应该没问题）。 T 介于 0.3 和 3.2 之间。我很确定问题在于b=T 存在一个点，因为我不断收到错误ValueError: Residuals are not finite in the initial point。但我不知道如何解决这个问题。

【问题讨论】：

永远不要仅仅假设数据是好的。即使在另一个地方没问题，也可能会导致这里出现问题。最好直接在示例中创建人工数据。这使得其他人可以验证问题（此外，找出导致问题的特定数据组合可能会引导您找到实际的解决方案）。
注意事项：如果T 介于0.3 和3.2 之间，如果b 介于1.8 和2.8 之间，您认为log(b-T) 的结果是什么？（此外，我不确定初始值是否全为 1，无论边界如何。）
@kazemakase 这很重要。无论如何，设置p0 可能是个好主意。此外，如果min(b) < max(T) 会遇到问题。所以应该检查一下。最后一个问题：为什么不拟合指数版本呢？这样可以避免log( negativNumber ) 问题。
@kazemakase 感谢您的意见。我对在这里发布问题有点陌生，所以我不太熟悉礼仪，所以我会记住你所说的关于数据的内容。
@IlinKaragjozov 这部分是礼仪，但主要是为了您自己的利益。你越容易让别人重现和解决你的问题，你就越有可能得到答案。如果可以简单地复制、粘贴和运行代码，人们就更有可能详细研究问题。

标签： python-3.x numpy scipy curve-fitting

【解决方案1】：

残差在初始点不是有限的

表示初始点不好，其中一些对数是无限的或未定义的。你需要一个更好的初始点。

根据模型的性质，b 必须大于 T 中的任何点。您目前拥有的b 上的边界不能保证这一点。把它们收紧。

当您不提供p0 参数时，SciPy 将在提供的范围内进行猜测。因此，如果边界保证有限性，则不会发生错误。不过，通常最好自己开p0，因为您对问题的先验理解比 SciPy 更好。

调整边界的工作示例：

popt, pcov=curve_fit(logfunc, np.linspace(0.3, 3.2, 6), [8, 7, 6, 5, 4, 3], bounds=([0.1, 3.2, 0.1], [1.0, 3.6, 1.0]))

【讨论】：

【解决方案2】：

您可能会发现 lmfit 包 (http://lmfit.github.io/lmfit-py/) 对于此类问题很有用。与scipy.optimize 包或curve_fit() 函数相比，这为曲线拟合问题提供了更高级别的方法以及对参数和模型的更好抽象。

对于这里的问题，lmfit的两个重要特点是

设置变量界限的能力。 curve_fit() 也可以做到这一点，但只能使用最小/最大边界的有序列表。使用lmfit，边界属于 Parameter 对象。
有一种方法可以显式设置处理 NaN 值的策略，这肯定会给您带来问题。

使用 lmfit，您的脚本将大致编写为

import numpy as np
import matplotlib.pyplot as plt
from lmfit import Model

def logfunc(T, a, b, c):
    return (a*np.log(b-T))+c

log_model = Model(logfunc, nan_policy='raise')  # raise error on NaNs
params = log_model.make_params(a=0.5, b=2.0, c=0.5) # initial values
params['b'].min = 1.8  # set min/max values
params['b'].max = 2.6 
params['c'].min = 0.1  # and so forth 

result = log_model.fit(np.log(Energy), params, T=T)

print(result.fit_report())

plt.plot(T, Energy, 'bo', label='data')
plt.plot(T, np.exp(result.best_fit), 'r--', label='fit')
plt.legend()
plt.xlabel('T')
plt.ylabel('Energy')
plt.gca().set_yscale('log', basey=10)
plt.show()

这比您的起始脚本稍微冗长，因为它提供了一个带标签的图，并且因为使用参数对象而不是标量提供了更大的灵活性和清晰度。

根据您的需要，您可以考虑将nan_policy 设置为“省略”，这将在出现NaN 时忽略它们——这不是一个好主意，但有时有助于您开始查找log(b-T) 的有效位置。你也可以改变你的模型函数来做类似的事情

def logfunc(T, a, b, c):
    arg = b - T
    arg[np.where(arg < 1.e-16)] = 1.e-16
    return a*np.log(arg) + c

明确防止出现 NaN 的一个明显原因。

【讨论】：