使用Minimumsq在python中同时进行数据拟合 [英] Simultaneous data fitting in python with leastsq

查看:118
本文介绍了使用Minimumsq在python中同时进行数据拟合的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我很久没有编程,也不擅长编程,但这是我正在努力的重要任务。我正在尝试拟合两组数据(x –时间,y1和y2 –应从文本文件中读取的不同值列)。对于每个数据集(y1和y2),我都有一个适合它们的函数。在这两个函数中,我都有几个要拟合的参数。
对于某些时间值, y的数据不存在,因此任务是当 y缺失时以某种方式对其进行编程,并在没有这些值的情况下拟合数据。参数应同时适合两个功能,因此应同时适合。那就是我无法解决的问题。

I didn't program for a long time and never was good at it, but it is kind of important task I am struggling with. I am trying to fit two sets of data (x – time, y1 and y2 – different columns of values which should be read from text file). For each dataset (y1 and y2) I have a function which should fit them. Inside both of the functions I have several parameters to be fitted. For some time values the data for "y" is absent, so the task is to program it somehow when "y" is missing and fit the data without these values. Parameters should fit both functions, so it should be simultaneous fitting. And that's the problem I cannot solve.

要在我的情况下定义函数,必须使用Inverse Laplace Transform,所以我使用了,对此没有任何疑问。

To define the function in my case it is necessary to use Inverse Laplace Transform, so I used and there are not questions about that.

import numpy as np
import pylab as plt
from scipy.optimize import leastsq
from cmath import *

# Here is the Laplace functions 

def Fp(s, td, m0, kon, koff):
   gs=s+kon-kon*koff/(s+koff)
   sr=np.sqrt(gs*td)
   return (m0/(s*s))*sr/sinh(sr)

def Fd(s, td, m0, kon, koff):
   gs=s+kon-kon*koff/(s+koff)
   sr=np.sqrt(gs*td)
   fu=koff/(koff+kon)
   fs=fu+koff*(1-fu)/(s+koff)
   return (m0/s)*fs*2*tanh(0.5*sr)/sr

# Define the trig functions cot(phi) and csc(phi)

def cot(phi):
   return 1.0/tan(phi)
def csc(phi):
   return 1.0/sin(phi)

# inverse Laplace transform for two functions which are going to be fitted next

def Qpt(t, td, m0, kon, koff): 
    shift = 0.1;
    ans = 0.0;
    N=30
    h = 2*pi/N;
    c1 = 0.5017
    c2 = 0.6407
    c3 = 0.6122
    c4 = 0. + 0.2645j

    for k in range(0,N):
         theta = -pi + (k+1./2)*h;
         z = shift + N/t*(c1*theta*cot(c2*theta) - c3 + c4*theta); 
         dz = N/t*(-c1*c2*theta*(csc(c2*theta)**2)+c1*cot(c2*theta)+c4);
         ans += exp(z*t)*Fp(z, td, m0, kon, koff)*dz;
    return ((h/(2j*pi))*ans).real

def Qdt(t,td, m0, kon, koff): 
    shift = 0.1;
    ans = 0.0;
    N=30
    h = 2*pi/N;
    c1 = 0.5017
    c2 = 0.6407
    c3 = 0.6122
    c4 = 0. + 0.2645j
    for k in range(0,N):
        theta = -pi + (k+1./2)*h;
        z = shift + N/t*(c1*theta*cot(c2*theta) - c3 + c4*theta); 
        dz = N/t*(-c1*c2*theta*(csc(c2*theta)**2)+c1*cot(c2*theta)+c4);
        ans += exp(z*t)*Fd(z, td, m0, kon, koff)*dz;
    return ((h/(2j*pi))*ans).real

#now we have Qp and Qd as theoretical functions

我编译并询问程序几个值,Qp和Qd正确定义。上面有关零件的唯一问题是:是否可以通过某种方式同时定义两个函数而无需执行两次转换?

I compiled that and asked a program several values, Qp and Qd are defined correctly. The only question about part above: is it possible somehow to define both functions at the same time without doing transform twice?

然后,我在零件中添加了同时拟合该函数和它给了我一个错误:

Then I added part with simultaneous fitting of this functions and it gave me a mistake:

TypeError: only length-1 arrays can be converted to Python scalars 

所以这是我适合的部分:

So this is my fitting part:

# FITTING PART

def residuals(pars, t, qd, qp):
    td = np.array(pars[0])
    m0 = np.array(pars[1])
    kon = np.array(pars[2])
    koff = np.array(pars[3])
    diff1 = Qdt(t,td, m0, kon, koff) - qd
    diff2 = Qpt(t,td, m0, kon, koff) - qp
    return np.concatenate((diff1[np.where(qd!=-1)], diff2[np.where(qp!=-1)]))

# for both functions with all the values
t = np.array([0.5, 2, 5, 10, 15, 20, 30, 40, 60, 90, 120, 180])
qd = np.array([0.145043746,0.273566338,0.437829373,0.637962531,-1,0.898107567,-1,1.186340492,1.359184345,-1,1.480552058,1.548143954])
qp = np.array([-1,-1,0.002701867,0.006485195,0.014034067,-1,0.06650739,-1,0.309055933,0.645945584,1.000811933,-1])


# initial values 
par_init = np.array([1, 1, 1, 1])

best, cov, info, mesg, ier = leastsq(residuals, par_init, args=(t, qd, qp), full_output=True)

print(" best-fit parameters: ", best)


#for each function separately to plot them and fitted functions as well
xqd= [0.5, 2, 5, 10, 20, 40, 60, 120, 180]
xqp= [5, 10, 15, 30, 60, 90, 120]
yqd= [0.145043746, 0.273566338, 0.437829373, 0.637962531, 0.898107567, 1.186340492, 1.359184345, 1.480552058, 1.548143954]
yqp= [0.002701867, 0.006485195, 0.014034067, 0.06650739, 0.309055933, 0.645945584, 1.000811933]

tt=np.linspace(0,185,100)
qd_fit=Qdt(tt,best[0], best[1], best[2], best[3])
qp_fit=Qdp(tt,best[0], best[1], best[2], best[3])

plt.plot(xqd,yqd,'bD:',xqp,yqp,'r^:', tt,qd_fit,'b',tt,qp_fit,'r')

plt.grid()
plt.show()

我将不胜感激!我迫切需要解决这个错误!

I would appreciate any help and advices! I desperately need to get it of this mistake!

预先感谢!

推荐答案

要同时使多个模型函数适应不同的数据集,您需要让残差函数执行以下操作:

to fit multiple model functions to different data sets simultaneously, you need to have your residuals function do the following:


  1. 具有全部将变量参数放入第1个数组中-最小sq()的第二个参数

  2. 接收所有要拟合的数据(此处,x,y1和y2)-必须与Minimumsq()的 args参数

  3. 解压缩当前参数值数组(残差的第一个参数)

  4. 计算每个子残差具有适当参数值的数据集

  5. 将结果连接到单个1d数组以获取返回值

  1. have all the variable parameters into 1 nd-array -- the second argument to leastsq()
  2. receive all the data to be fitted (here, x, y1, and y2) -- must match the "args" argument of leastsq()
  3. unpack the "current parameter value" array (residuals 1st argument)
  4. calculate the sub-residual for each dataset with the values of the appropriate parameters
  5. concatenate the results to a single 1d-array for the return value

一个简单的残差可能看起来像这样:

A simple residual might look like this:

import numpy as np
from scipy.optimize import leastsq

def residual_two_functions(pars, x, y1, y2):
    off1   = pars[0]
    slope1 = pars[1]
    off2   = pars[2]
    slope2 = pars[3]
    diff1 = y1 - (off1 + slope1 * x)
    diff2 = y2 - (off2 + slope2 * x)
    return np.concatenate((diff1, diff2))

# create two tests data sets
NPTS = 201
x  = np.linspace(0, 10, NPTS)
y1 = -0.7 +  1.7*x + np.random.normal(scale=0.01, size=NPTS)
y2 =  5.2 + 50.1*x + np.random.normal(scale=0.02, size=NPTS)

# initial values 
par_init = np.array([0, 1, 10, 100])

best, cov, info, message, ier = leastsq(residual_two_functions,
                                        par_init, args=(x, y1, y2),
                                        full_output=True)

print(" Best-Fit Parameters: ",  best)

这篇关于使用Minimumsq在python中同时进行数据拟合的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

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