python-3.x 多指数衰减的lmfit模型

6qfn3psc  于 2023-11-20  发布在  Python
关注(0)|答案(2)|浏览(162)

我尝试使用lmfit包对我的数据进行拟合。但是我找不到任何内置的多指数衰减模型。我尝试创建自己的函数,然后拟合它。我的代码如下:

import os
import time
import numpy as np
import pandas as pd
import lmfit
from lmfit.models import ExponentialModel, LinearModel
from lmfit import Model, Parameter, report_fit

def MultiExpDecay(tiempo,C1,tau1,C2,tau2):
    return C1*np.exp(-tiempo/tau1)+C2*np.exp(-tiempo/tau2)

def MultiExpDecay_fit():
    C1s = []
    C1s_error = []
    C2s = []
    C2s_error = []
    tau1s = []
    tau1s_error = []
    tau2s = []
    tau2s_error = []
    Fit_MultiExpDecays = []
    model = Model(MultiExpDecay, independent_vars=['tiempo'])
    for c in range(len(V_APD.columns)):
        xdat = tiempo.iloc[:, c]
        ydat = V_APD.iloc[:, c]
        pars = model.guess(ydat, x=xdat)
        fit = model.fit(ydat, pars, x=xdat)
        fit_values = model.eval(pars, x=xdat)
        Fit_MultiExpDecays.append(fit_values)
        for key in fit.params:
            if key == 'C1':
                #print(key, "=", out.params[key].value, "+/-", out.params[key].stderr)
                C1s.append(fit.params[key].value)
                C1s_error.append(fit.params[key].stderr)
            elif key == 'C2':
                C2s.append(fit.params[key].value)
                C2s_error.append(fit.params[key].stderr)
            elif key == 'tau1':
                tau1s.append(fit.params[key].value)
                tau1s_error.append(fit.params[key].stderr)
            elif key == 'tau2':
                tau2s.append(fit.params[key].value)
                tau2s_error.append(fit.params[key].stderr)
    Fit_MultiExpDecays = np.transpose(pd.DataFrame(Fit_MultiExpDecays, index=labels))
    C1 = np.transpose(pd.DataFrame(C1s, index = labels, columns = ['C1']))
    C2 = np.transpose(pd.DataFrame(C2s, index=labels, columns=['C2']))
    tau1 = np.transpose(pd.DataFrame(tau1s, index=labels, columns=['tau1']))
    tau2 = np.transpose(pd.DataFrame(tau2s, index=labels, columns=['tau2']))
    C1_error = np.transpose(pd.DataFrame(C1s_error, index = labels, columns=['C1 error']))
    C2_error = np.transpose(pd.DataFrame(C2s_error, index=labels, columns=['C2 error']))
    tau1_error = np.transpose(pd.DataFrame(tau1s_error, index=labels, columns=['tau1 error']))
    tau2_error = np.transpose(pd.DataFrame(tau2s_error, index=labels, columns=['tau2 error']))
    C1 = pd.concat([C1, C1_error])
    C2 = pd.concat([C2, C2_error])
    tau1 = pd.concat([tau1, tau1_error])
    tau2 = pd.concat([tau2, tau2_error])
    return C1, C2, tau1, tau2, Fit_MultiExpDecays

C1, C2, tau1, tau2, Fit_MultiExpDecays = MultiExpDEcay_fit():

字符串
出现错误,但无法识别问题。

File "C:\ProgramData\Anaconda3\Lib\site-packages\lmfit\model.py", line 737, in guess
    raise NotImplementedError(msg)
NotImplementedError: guess() not implemented for Model

thigvfpy

thigvfpy1#

也可以使用两个lmfit.ExponentialModel

import numpy as np
from lmfit.models import ExponentialModel

def dblexp( x, c1, l1, c2, l2 ):
    return c1 * np.exp( -x / l1 ) + c2 * np.exp( -x / l2 )
    
xl = np.linspace(0, 10, 101)
yl = dblexp(xl, 32.44, 0.81, 11.51, 9.22 ) + np.random.normal(0, 0.2, xl.size)

mymodel = ExponentialModel(prefix='e1_') + ExponentialModel(prefix='e2_') 

params = mymodel.make_params(e1_amplitude=10, e1_decay=10,
                             e2_amplitude=25, e2_decay=0.50)

result = mymodel.fit(yl, params, x=xl)
print( result.fit_report() )

字符串
对于实际问题:正如异常所说,Model不知道如何为任意的用户定义的模型函数猜测合理的参数值-Model.guess()方法在默认情况下引发此异常,并且必须由每个子类覆盖。lmfit为提供的子模型(包括ExponentialModel)提供此方法,但这些方法都是在理解这些特定模型的情况下编码的。
最后:拟合双指数衰减是出了名的容易出错,并且对于最小二乘方法(或至少Levenberg-Marquardt实现)来说是困难的。对于更现实的情况,我建议用双指数衰减的对数来拟合数据的对数。

mnemlml8

mnemlml82#

lmfit的设计使你不必做你在这里做的所有编程。它应该看起来像:

import matplotlib.pyplot as plt
from numpy import linspace
from numpy import fromiter
from numpy import exp
from numpy.random import normal
from lmfit import Model

def dblexp( x, c1, l1, c2, l2 ):
    return c1 * exp( -x / l1 ) + c2 * exp( -x / l2 )
    
xl = linspace(0, 10, 101)
yl = dblexp( xl, 32.44, 0.81, 11.51, 9.22 ) + normal(0, 0.2, xl.size)

mymodel = Model( dblexp ) 

### need some good guesses to start with
params = mymodel.make_params(c1=30, l1=1, c2=5, l2 =10)

result = mymodel.fit(yl, params, x=xl)
print( result.fit_report() )

fig = plt.figure()
ax = fig.add_subplot( 1, 1, 1 )
ax.scatter( xl, yl )
ax.plot( xl, result.best_fit, 'r-' )
plt.show()

字符串
提供:

[[Model]]
    Model(dblexp)
[[Fit Statistics]]
    # fitting method   = leastsq
    # function evals   = 26
    # data points      = 101
    # variables        = 4
    chi-square         = 4.77426620
    reduced chi-square = 0.04921924
    Akaike info crit   = -300.239903
    Bayesian info crit = -289.779421
[[Variables]]
    c1:  32.5481660 +/- 0.18540661 (0.57%) (init = 30)
    l1:  0.79183180 +/- 0.00931390 (1.18%) (init = 1)
    c2:  11.6241194 +/- 0.16669056 (1.43%) (init = 5)
    l2:  9.11790608 +/- 0.19623957 (2.15%) (init = 10)
[[Correlations]] (unreported correlations are < 0.100)
    C(c2, l2) = -0.949
    C(l1, c2) = -0.822
    C(l1, l2) =  0.733
    C(c1, c2) = -0.652
    C(c1, l2) =  0.645
    C(c1, l1) =  0.251



的数据

相关问题