1

我正在尝试重现已知的拟合结果(在期刊论文中报告):将幂律模型应用于数据。从下面的图 A 中可以看出,我能够通过使用已知的最佳拟合参数来重现结果。

< Plot-A:来自文献的已知结果 > 两个轴:对数刻度

但是,我无法自己重新推导最佳拟合参数。

< Plot-B:curve_fit 和 lmfit 的不正确拟合 > 错误的

案例A返回,

OptimizeWarning: Covariance of the parameters could not be estimated(如果我省略了几个初始数据点,拟合会返回一些不错的结果,但仍与已知的最佳拟合结果不同)。 编辑:现在,我这次只是发现了新的附加错误消息.. : (1) RuntimeWarning: overflow encountered in power (2)RuntimeWarning: invalid value encountered in power

case-B(初始猜测更接近最佳拟合参数)返回,

RuntimeError: Optimal parameters not found: Number of calls to function has reached maxfev = 5000. 如果我将 maxfev 设置得更高以考虑此错误消息,则拟合有效但返回不正确的结果(与最佳拟合结果相比非常错误的拟合)。

4

1 回答 1

2

正如我评论的那样:

由于您将数据绘制在对数图上,您是否也适合 log(y) 和 log(x)?由于您的 y 数据变化 5 或 6 个数量级,如果您不适合对数空间,则只有具有最高 y 值的 3 或 4 个数据点才重要。

显然,这暗示有点太微妙了。所以我会更直接:如果您在日志空间中绘图,则适合日志空间。

而且,您的模型很容易从和 NaN 生成复数negative**fraction,这无疑会导致您的所有拟合出现问题。总是打印出最佳拟合参数和协方差矩阵。

因此,您可能需要对参数施加限制(当然,我不知道您的模型是否正确,或者实际上是您认为“正确答案”使用的)。也许从这样的事情开始:

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

# the data can be found at the end of this post.
xx, yy = np.genfromtxt('the_data', unpack=True)

# the BPL function
def bendp(x, A, x_bend, allo, alhi, c):
    numer = A * x**-allo
    denom = 1 + (x/x_bend)**(alhi-allo)
    out = (numer/denom) + c
    return  np.log(out)       ## <- TAKE THE LOG


mod = Model(bendp)
para = mod.make_params(A = 0.01, x_bend=1e-2, allo=1., alhi=2., c=1e-2)

# Note: your model is very sensitive # make sure A, x_bend, alhi, and c cannot be negative
para['A'].min = 0
para['x_bend'].min = 0
para['alhi'].min = 0
para['alhi'].max = 3
para['c'].min = 0


result = mod.fit(np.log(yy), para, x=xx) ## <- FIT THE LOG

print(result.fit_report())

plt.loglog(xx, yy, linewidth=0, marker='o', markersize=6)
plt.loglog(xx, np.exp(result.best_fit), 'r', linewidth=5)
plt.show()

希望有帮助...

于 2020-01-15T04:47:04.830 回答