如何使用python拟合多指数曲线 [英] How to fit multiple exponential curves using python

查看:34
本文介绍了如何使用python拟合多指数曲线的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

对于这里的图像curve_fit for as single exponential curve所示的单个指数曲线,我可以使用scipy.Optimize.curveFit来拟合数据。然而,我不确定如何实现对这里所示的由多个指数曲线组成的相似数据集的拟合double exponential curves。 我使用以下方法实现了对单曲线的拟合:

def exp_decay(x,a,r):
    return a * ((1-r)**x) 

x = np.linspace(0,50,50)
y = exp_decay(x, 400, 0.06)

y1 = exp_decay(x, 550, 0.06)      # this is to be used to append to y to generate two curves

pars, cov = curve_fit(exp_decay, x, y, p0=[0,0])
plt.scatter(x,y)
plt.plot(x, exp_decay(x, *pars), 'r-')     #this realizes the fit for a single curve

yx = np.append(y,y1)   #this realizes two exponential curves (as shown above - double exponential curves) for which I don't need to fit a model to

有人可以帮助描述一下如何为包含两条曲线的数据集实现这一点吗?我的实际数据集由多条指数曲线组成,但我认为如果我能够实现对两条曲线的拟合,我也许能够为我的数据集复制相同的曲线。这不能使用Scipy的curve_fit来完成;任何可以工作的实现都可以。

请帮帮忙!

推荐答案

您的问题可以通过使用一阶导数估计等简单标准拆分您的数据集轻松解决,然后我们可以对每个子数据集应用简单的曲线拟合过程。

试验数据集

首先,让我们导入一些包并创建一个包含三条曲线的合成数据集来表示您的问题。

我们使用一个双参数指数模型,因为时间原点偏移将由分裂方法处理。我们还会添加噪声,因为实际数据中始终存在噪声:

import numpy as np
import pandas as pd
from scipy import optimize
import matplotlib.pyplot as plt

def func(x, a, b):
    return a*np.exp(b*x)

N = 1001
n1 = N//3
n2 = 2*n1

t = np.linspace(0, 10, N)

x0 = func(t[:n1], 1, -0.2)
x1 = func(t[n1:n2]-t[n1], 5, -0.4)
x2 = func(t[n2:]-t[n2], 2, -1.2)

x = np.hstack([x0, x1, x2])
xr = x + 0.025*np.random.randn(x.size)

以图形方式呈现如下:

数据集拆分

我们可以将数据集分成三个子数据集,使用一个简单的准则作为一阶导数估计,使用一阶差分对其进行评估。目标是检测曲线何时急剧上升或下降(数据集应该拆分的位置)。一阶导数估计如下):

dxrdt = np.abs(np.diff(xr)/np.diff(t))

标准需要额外的参数(阈值),必须根据您的信号规范进行调整。条件相当于:

xcrit = 20
q = np.where(dxrdt > xcrit) # (array([332, 665], dtype=int64),)

和拆分指数为:

idx = [0] + list(q[0]+1) + [t.size] # [0, 333, 666, 1001]

标准阈值主要受数据上噪声的性质和强度以及两条曲线之间的间隙大小的影响。此方法的使用取决于在存在噪声的情况下检测曲线间隙的能力。当噪声功率与我们要检测的间隙的大小相同时,它就会破裂。如果噪声严重拖尾(很少是强异常值),您还可以观察到错误的分裂指数。

在此MCVE中,我们已将阈值设置为20 [Signal Units/Time Units]

此手工创建的标准的替代方案是将标识委托给scipy的优秀find_peaks方法。但它无法避免根据您的信号规格调整检测的要求。

拟合原点偏移的数据集

现在我们可以对每个子数据集(原点偏移时间)进行曲线拟合,收集参数和统计数据,并绘制结果:

trials = []
fig, axe = plt.subplots()
for k, (i, j) in enumerate(zip(idx[:-1], idx[1:])):
    p, s = optimize.curve_fit(func, t[i:j]-t[i], xr[i:j])
    axe.plot(t[i:j], xr[i:j], '.', label="Data #{}".format(k+1))
    axe.plot(t[i:j], func(t[i:j]-t[i], *p), label="Data Fit #{}".format(k+1))
    trials.append({"n0": i, "n1": j, "t0": t[i], "a": p[0], "b": p[1],
                   "s_a": s[0,0], "s_b": s[1,1], "s_ab": s[0,1]})
axe.set_title("Curve Fits")
axe.set_xlabel("Time, $t$")
axe.set_ylabel("Signal Estimate, $hat{g}(t)$")
axe.legend()
axe.grid()
df = pd.DataFrame(trials)

它返回以下拟合结果:

    n0    n1    t0         a         b       s_a           s_b      s_ab
0    0   333  0.00  0.998032 -0.199102  0.000011  4.199937e-06 -0.000005
1  333   666  3.33  5.001710 -0.399537  0.000013  3.072542e-07 -0.000002
2  666  1001  6.66  2.002495 -1.203943  0.000030  2.256274e-05 -0.000018

符合我们的原始参数(参见试验数据集部分)。

我们可以用图形方式检查拟合优度:

这篇关于如何使用python拟合多指数曲线的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆