2016-04-21 42 views
2

我对使用python进行数据分析以及将lmfit用于非线性方程拟合有点新。我试图建立一个复杂的半解析函数,描述一维通道中污染物羽流随时间推移的运动。我正在使用lmfit v0.9.3。我已经成功地通过了lmfit教程中的几个示例,但似乎无法让我自己的模型起作用。下面的脚本,作品达tsm_mod.fit()呼叫,但随后返回错误:使用lmfit的错误v0_9_3

File "C:\Anaconda\lib\site-packages\lmfit\model.py", line 501, in fit for p in params.values()])

TypeError: 'numpy.ndarray' object is not callable

的代码如下:

import numpy as np 
import scipy as sp 
import pandas as pd 
from lmfit import Model 

desmedt = pd.read_table('Directory\desmedt_test.txt',sep='\t') 

x = desmedt['Times'] 
y = desmedt['Conc'] 

def tsm_intfunc(t,x,tau,u,k,alpha,beta,mass,ac): 
    return((mass/(2*ac*(t*np.pi*k)**(1/2)))*np.exp(-((x-u*t)**2)/(4*k*t))*np.exp(-alpha*tau-alpha*(t-tau)/beta) 
     *np.sqrt(beta*tau/(t-tau))*sp.special.iv(2*np.sqrt((alpha**2)*tau*(t-tau)/beta),1)) 

def tsm_desmedt(t,x,u,k,alpha,beta,mass,ac,nsteps): 
    dtau = t/nsteps 
    cxt = (mass/(2*ac*np.sqrt(t*np.pi*k)))*np.exp(-((x-u*t)**2)/(4*k*t))*np.exp(-alpha*t) 
    cxv = tsm_intfunc(t,x,0.00000001,u,k,alpha,beta,mass,ac)/2 
    i = 1 

    while (i<nsteps): 
     cxv = cxv+tsm_intfunc(t,x,dtau*i,u,k,alpha,beta,mass,ac)/2 
     i = i+1 

    return cxt+(alpha/beta)*cxv*dtau 

tsm_mod = Model(tsm_desmedt) 
tsm_mod.set_param_hint('ac',value=18.2,vary=False) 
tsm_mod.set_param_hint('alpha',value=1e-4) 
tsm_mod.set_param_hint('beta',value=1e-1) 
tsm_mod.set_param_hint('k',value=3) 
tsm_mod.set_param_hint('mass',value=157100,vary=False) 
tsm_mod.set_param_hint('nsteps',value=100,vary=False) 
tsm_mod.set_param_hint('u',value=0.4) 
tsm_mod.set_param_hint('x',value=4604,vary=False) 

tsm_pars = tsm_mod.make_params() 
tsm_fit = tsm_mod.fit(y,x,tsm_pars) 

难道这是在lmfit的错误吗?或者,您认为我使用lmfit设置问题的方式有错误吗?

编辑: 在装配中使用的数据在下面给出:

时报

7787.628 8330.04 8640 8756.244 8988.696 9143.676 9337.392 9492.372 9724.86 9918.576 10034.784 10228.536 10383.516 10577.232 10770.948 11003.4 11119.644 11313.36 11468.34 11700.792 11855.772 12010.752 12204.468 12359.448 12630.672 12824.388 13173.084 13483.044 13793.004 14412.924 14955.336 15575.256 16195.14 17357.472

0.00944669 0.0850202 0.236167 0.576248 1.00135 2.01215 2.84345 3.51417 4.53441 5.21457 5.59244 5.74359 5.88529 6.0081 5.75304 5.61134 5.20513 4.95007 4.41161 3.74089 3.46694 3.07962 2.80567 2.41835 2.1444 1.74764 1.47368 1.20918 0.935223 0.661269 0.406208 0.132254 0.11336 0.151147

+0

你能提供你正在使用的数据还是至少它的一个子集? – Cleb

+0

根据[docs online](https://lmfit.github.io/lmfit-py/model.html#model.Model),'Model.fit'预计会传递参数,参数,权重,.. 。你已经给它y,x,tsm \ _pars,这看起来不太合适。特别是,它会将你的'x'视为它的'params',它预计它是'Parameters'对象。 –

+0

只要我在网上看到的东西 - 我从来没有用过py-lmfit - 我想也许你想要像'tsm_fit = tsm_mod.fit(y,tsm_pars,x = x)'这样的东西。 –

回答

0

你的例子稍微在这混乱你的模型函数包含一个独立变量t和一个变量x,并在您的主程序中调用时间x。无论如何,你想要的电话是

tsm_fit = tsm_mod.fit(y, tsm_pars, t=x) 

该拟合似乎不太好,但这似乎是一个单独的问题。

+0

是的,谢谢你的意见。我意识到了代码中的错误,并能够成功实现它(请参阅上面的注释)。然而,该模型具有参数可识别性问题,因此如果初始参数猜测远离其最佳值,则该拟合较差。有趣的是,我也在R中实现了这个算法,并且能够更好地适应(使用LMFit包),使参数值更接近他们应该在的时间比使用Python更少。最有趣的是Excel解算器给出了参数值的最佳估计,但不能报告标准错误。 –