Fit Multiple Parametric Curves With Scipy

I have a set (at least 3) of curves (xy-data). For each curve the parameters E and T are constant but different. I'm searching the coefficients a,n and m for the best fit over all

Solution 1:

thanks Evert for the reply.

Exactly what I needed to know!!

I have simplyfied the function as far as possible, as you suggested. However, the task was to find ONE set of A,m,n to fit all curves. So my code look like this:

import numpy
import math
from scipy.optimize import leastsq

#+++++++++++++++++++++++++++++++++++++++++++++deffit(x,T,A,n,m):
  return A/(n+1.0)*math.pow(T,(n+1.0))*numpy.power(x,m)
#+++++++++++++++++++++++++++++++++++++++++++++defleastsq_func(params, *args):
  cc=args[0]   #number of curves
  incs=args[1] #number of points 
  x=args[2]
  y=args[3]
  T=args[4:]

  A=params[0]
  n=params[1]
  m=params[2]

  yfit=numpy.empty(x.shape)
  for i inrange(cc):
    v=i*incs
    b=(i+1)*incs
    if b<cc:
     yfit[v:b]=fit(x[v:b],T[i],A,n,m)
    else:
     yfit[v:]=fit(x[v:],T[i],A,n,m)

  return y-yfit
#+++++++++++++++++++++++++++++++++++++++++++++
Ts  =[10,100,1000,10000]    #4 T-values for 4 curves
incs=10#10 datapoints in each curve
x=["measured data"]   #all 40 x-values
y=["measrued data"]   #all 40 y-values
x=numpy.array(x)   
y=numpy.array(y)   

params0=[0.001,1.01,-0.8]   #parameter guess

args=[len(Ts),incs,x,y]
for c in Ts:
  args.append(c) 
args=tuple(args)   #doesn't work if args is a list!!

result=leastsq(leastsq_func, params0, args=args)

Works like clockwork.

At first I put the Ts in the params0 list and they were modified during iteration leading to nonsense results. Obvious, if you think about it. Afterwards ;-)

So, Vielen Dank! J.

Solution 2:

Thank you guys, I found this very useful. In case someone wants a general solution to this problem, I wrote one that is heavily inspired by the snippets above:

import numpy as np
from scipy.optimize import leastsq

defmultiple_reg(x, y, f, const, params0, **kwargs):
    """Do same non-linear regression on multiple curves
    """defleastsq_func(params, *args):
        x, y = args[:2]
        const = args[2:]
        yfit = []
        for i inrange(len(x)):
            yfit = np.append(yfit, f(x[i],*const[i],*params))
        return y-yfit

    # turn const into 2d-array if 1d is given
    const = np.asarray(const)
    iflen(const.shape) < 2:
        const = np.atleast_2d(const).T

    # ensure that y is flat and x is nestedifhasattr(y[0], "__len__"):
        y = [item for sublist in y for item in sublist]
    ifnothasattr(x[0], "__len__"):
        x = np.tile(x, (len(const), 1))
    x_ = [item for sublist in x for item in sublist]
    assertlen(x_) == len(y)

    # collect all arguments in a tuple
    y = np.asarray(y)
    args=[x,y] + list(const)
    args=tuple(args)   #doesn't work if args is a list!!return leastsq(leastsq_func, params0, args=args, **kwargs)

This function accepts xs and ys of various length, as they are stored in nested lists rather than numpy ndarrays. For the particular case presented in this thread, the function can be used like this:

def fit(x,T,A,n,m):
    return A/(n+1.0)*np.power(T,(n+1.0))*np.power(x,m)

# prepare dataset with some noise
params0 = [0.001, 1.01, -0.8]
Ts = [10, 50]

x = np.linspace(10, 100, 100)
y = np.empty((len(Ts), len(x)))
for i in range(len(Ts)):
    y[i] = fit(x, Ts[i], *params) + np.random.uniform(0, 0.01, size=len(x))

# do regression
opt_params, _ = multiple_reg(x, y, fit, Ts, params0)

Verify regression by plotting the data and regression lines

import matplotlib.pyplot as plt
for i in range(len(Ts)):
    plt.scatter(x, y[i], label=f"T={Ts[i]}")
    plt.plot(x, fit(x, Ts[i], *opt_params), '--k')
plt.legend(loc='best')
plt.show()