【问题标题】:Least squares function and 4 parameter logistics function not working最小二乘函数和 4 参数物流函数不起作用
【发布时间】:2017-10-13 03:13:04
【问题描述】:

python 相对较新,主要用于绘图。我目前正在尝试使用 4 参数逻辑 (4PL) 方程和来自 scipy 的曲线拟合来确定最佳拟合线。有一两个网站展示了 4PL 的工作原理,但无法让它们为我的数据工作。示例,但类似的 4PL 数据如下:

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
import scipy.optimize as optimization

xdata = [2.3, 2.3, 2, 2, 1.7, 1.7, 1, 1, 0.000001, 0.000001, -1, -1]
ydata = [0.32, 0.3, 0.55, 0.60, 0.88, 0.92, 1.27, 1.21, 1.15, 1.12, 1.1, 1.1]
def fourPL(x, A, B, C, D):
    return ((A-D)/(1.0+((x/C)**(B))) + D)
guess = [0, -0.5, 0.5, 1]
params, params_covariance = optimization.curve_fit(fourPL, xdata, ydata, 
guess)
params

给出警告(也是测试数据中的指数警告,但不是真实的):

OptimizeWarning: Covariance of the parameters could not be estimated
category=OptimizeWarning)

参数返回我最初的猜测。我尝试了各种初步猜测。

在绘图时绘制了最佳拟合线,但它不是曲线并且不低于 x = 0(我找不到底片会与 4PL 模型混淆的原因)。 4PL fit plotted

我不确定我是否对方程做错了,或者曲线拟合函数如何工作,或两者兼而有之。我有一个使用最小二乘法而不是曲线拟合的类似问题。我已经尝试了一系列基于适合等式的类似方程的变化,但已经卡住了一段时间,任何帮助我指向正确方向的帮助将不胜感激。

【问题讨论】:

  • 您好,欢迎您!请提供minimal, complete, and verifyable example(即带有玩具数据的完整代码)以增加获得有用答案的机会。人们只想简单地复制和粘贴您的代码以重现错误并摆弄可能的解决方案。也看看here。我昨天发布的这个答案可能会给你一些想法,尽管逻辑函数应该很适合curve_fit
  • 我从图中看到,对于图中的每个 X 值,您都有多个数据点。它可能会简化您对这些值进行平均的初始工作,从而在您开始时为您提供一个较小的数据集。一旦您使用较小的数据集获得了良好的工作解决方案,然后切换回使用较大的完整数据集。这可能是一种快速简便且易于执行的测试。
  • 感谢您的建议,编辑后有一个简单、可复制的示例。还玩弄了示例数据,似乎卡在了同一个地方。和以前一样。

标签: python scipy curve-fitting


【解决方案1】:

我很惊讶您没有收到任何警告或没有与我们分享。我无法通过科学的方式为您分析此任务,仅对技术内容进行一些说明:

观察

在运行您的代码时,您应该出现如下警告:

RuntimeWarning: invalid value encountered in power
  return ((A-D)/(1.0+((x/C)**(B))) + D)

不要忽视这一点!

调试

为你的函数fourPL添加一些打印,可能是你函数的所有不同组件,看看发生了什么。

例子:

def fourPL(x, A, B, C, D):
    print('1: ', (A-D))
    print('2: ', (x/C))
    print('3: ', (1.0+((x/C)**(B))))

    return ((A-D)/(1.0+((x/C)**(B))) + D)

...
params, params_covariance = optimization.curve_fit(fourPL, xdata, ydata, guess, maxfev=1)
# maxfev=1 -> let's just check 1 or few it's

输出:

1:  -1.0
2:  [  4.60000000e+00   4.60000000e+00   4.00000000e+00   4.00000000e+00
   3.40000000e+00   3.40000000e+00   2.00000000e+00   2.00000000e+00
   2.00000000e-06   2.00000000e-06  -2.00000000e+00  -2.00000000e+00]
RuntimeWarning: invalid value encountered in power
  print('3: ', (1.0+((x/C)**(B))))
3:  [   1.4662524     1.4662524     1.5           1.5           1.54232614
    1.54232614    1.70710678    1.70710678  708.10678119  708.10678119
           nan           nan]

这就够了。 nans 和 infs 很糟糕!

理论

现在是理论的时候了,我不会那样做。但通常你现在应该考虑一下基本理论以及为什么会出现这些问题。

关于假设,您是否遗漏了什么?

修理(不检查理论)

没有检查理论,只是查看在 30 秒内发现的一些 example:嗯,负 x 值有问题吗?

让我们移动 x(最小;这里硬编码 1):

xdata = np.array([2.3, 2.3, 2, 2, 1.7, 1.7, 1, 1, 0.000001, 0.000001, -1, -1]) + 1

完整代码:

import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
import scipy.optimize as optimization

xdata = np.array([2.3, 2.3, 2, 2, 1.7, 1.7, 1, 1, 0.000001, 0.000001, -1, -1]) + 1
ydata = np.array([0.32, 0.3, 0.55, 0.60, 0.88, 0.92, 1.27, 1.21, 1.15, 1.12, 1.1, 1.1])

def fourPL(x, A, B, C, D):
    return ((A-D)/(1.0+((x/C)**(B))) + D)

guess = [0, -0.5, 0.5, 1]
params, params_covariance = optimization.curve_fit(fourPL, xdata, ydata, guess)#, maxfev=1)

x_min, x_max = np.amin(xdata), np.amax(xdata)
xs = np.linspace(x_min, x_max, 1000)
plt.scatter(xdata, ydata)
plt.plot(xs, fourPL(xs, *params))
plt.show()

输出:

RuntimeWarning: divide by zero encountered in power
  return ((A-D)/(1.0+((x/C)**(B))) + D)

看起来不错,但现在是另一场理论课的时候了:我们的线性位移对我们的结果有什么影响?我又忽略了。

所以只有一个警告和漂亮的输出。

如果您想删除最后一条警告,请添加一些小 epsilon 以使 xdata 中没有 0:

xdata = np.array([2.3, 2.3, 2, 2, 1.7, 1.7, 1, 1, 0.000001, 0.000001, -1, -1]) + 1 + 1e-10

这将达到相同的效果,没有任何警告。

【讨论】:

  • 谢谢萨沙!效果很好。我想我可以只绘制我的数据而不首先转换为日志,这使所有数字保持正数和非零。
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2011-09-16
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多