scipy 的指数衰减只给出阶跃函数

exponential decay with scipy just gives step function

我正在尝试对一组数据进行指数拟合:

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


def func(x, a, b, c):
    return a * np.exp(x / -b) + c


epr_data = np.loadtxt('T2_text', skiprows=1)

time = epr_data[:, 1]
intensity = epr_data[:, 2]


optimizedParameters, pcov = opt.curve_fit(func, time, intensity)
print(optimizedParameters)
plt.plot(time, intensity, func(time, *optimizedParameters), label="fit")

plt.show()

但我只得到这个阶跃函数和这些参数:

[1.88476367e+05 1.00000000e+00 6.49563230e+03]

“适合”的情节

以及此错误消息:

OptimizeWarning: Covariance of the parameters could not be estimated
  warnings.warn('Covariance of the parameters could not be estimated'

编辑:

https://pastebin.com/GTTGf0ed

我想用强度绘制时间和第一行

您建议后的图表

编辑 2:

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


def func(x, a, b, c):
    return a * np.exp(x / -b) + c


epr_data = np.loadtxt('T2_text', skiprows=1)

time = epr_data[:, 1]
intensity = epr_data[:, 2]

c0 = np.mean(intensity[-10:])
a0 = intensity[0]-c0
th = time[np.searchsorted(-intensity+c0, -0.5*a0)]
b0 = th / np.log(2)

optimizedParameters, pcov = opt.curve_fit(func, time, intensity, p0=(a0, b0, c0))
print(optimizedParameters)
plt.plot(time, intensity, label='data')
plt.plot(time, func(time, *optimizedParameters), label="fit")
plt.legend()

plt.show()

首先,修复您的 plot 函数调用。要一次调用绘制两条曲线,您必须为每条曲线提供 xy 值:

plt.plot(time, intensity, time, func(time, *optimizedParameters), label="fit")

为了标注,调用plot两次可能更简单:

plt.plot(time, intensity, label="data")
plt.plot(time, func(time, *optimizedParameters), label="fit")

如果我们有您的数据,解决 curve_fit 生成的警告会更容易。然而,经验表明,curve_fit 使用的参数的默认初始猜测(全为 1)很可能只是一个非常糟糕的猜测。

尝试通过为其数值优化例程提供更好的起点来帮助 curve_fit。下面的代码显示了如何计算 abc 的粗略猜测。使用参数 p0=(a0, b0, c0).

将这些传递给 curve_fit
# Assume that the time series is long enough that the tail
# has flattened out to approximately random noise around c.
c0 = np.mean(intensity[-10:])

# This assumes time[0] is 0.
a0 = intensity[0] - c0

# Rough guess of the half-life.
th = time[np.searchsorted(-intensity+c0, -0.5*a0)]
b0 = th / np.log(2)