Search code examples
pythonlmfit

Using Model class in lmfit in combination with Parameter class


I'm trying to fit some data using lmfit. Originally, I used the Model class which worked well based on this example/tutorial: https://lmfit.github.io/lmfit-py/model.html But then I wanted to add some parameters constraints to the model, thus I looked at this tutorial: https://lmfit.github.io/lmfit-py/parameters.html

However, I've got some problems in combining the two classes such that they work nicely together. Either it complains about the fit function missing parameters or getting invalid parameters (that's the case in the example I'll post) or I get a model which doesn't actually take the parameters I specified. I can actually solve the problem using one of these different approaches: 1. Pass the paramaters using model.make_params(...), but I would like to split them up individually 2. I could use the Minimizer instead of Model, but I would like to understand why they are so differently implemented even though I would expect them to be very similar (except that they work on different type of input)

Any help / explanations would be really appreciated. :)

This code is based on the example on the Parameters tutorial page. What I did here is to modify the example such that the Model class is used instead of the Minimizer class, which in principle should work, but I'm somehow doing it the wrong way. As comparison, the original example is here (scroll to the bottom): https://lmfit.github.io/lmfit-py/parameters.html

# <examples/doc_parameters_basic.py>
import numpy as np

from lmfit import Model, Parameters

# create data to be fitted
x = np.linspace(0, 15, 301)
data = (5. * np.sin(2*x - 0.1) * np.exp(-x*x*0.025) +
        np.random.normal(size=len(x), scale=0.2))


# define objective function: returns the array to be minimized
def fcn2min(params, x):
    """Model a decaying sine wave and subtract data."""
    amp = params['amp']
    shift = params['shift']
    omega = params['omega']
    decay = params['decay']
    model = amp * np.sin(x*omega + shift) * np.exp(-x*x*decay)
    return model


# create a set of Parameters
params = Parameters()
params.add('amp', value=10, min=0)
params.add('decay', value=0.1)
params.add('shift', value=0.0, min=-np.pi/2., max=np.pi/2)
params.add('omega', value=3.0)

# do fit, here with leastsq model
fitmodel = Model(fcn2min)
result = fitmodel.fit(data, params, x=x)

# calculate final result
final = result.fit_report()

# try to plot results
try:
    import matplotlib.pyplot as plt
    plt.plot(x, data, 'k+')
    plt.plot(x, final, 'r')
    plt.show()
except ImportError:
    pass
# <end of examples/doc_parameters_basic.py>

But using it like this, I get an error

ValueError: Invalid independent variable name ('params') for function fcn2min

I tried: Specifying all parameters as function arguments, e.g.

def fcn2min(x, amp, shift, omega, decay):

But in this case I end up with the model / function parameters not being connected and I get a 'fit' which does nothing.

Now I tried stuff like specifying:

fitmodel = Model(fcn2min, independent_vars=['x'], param_names=params)

But in this case I get

Invalid parameter name ('amp') for function fcn2min

Also tried something like this:

params = fitmodel.make_params()
params.add('amp', value=10, min=0)
...

But in this case I also don't get parameters which are connected to the model, which can be seen from the output of

fitmodel.param_names

which returns an empty list.


Solution

  • You are confusing a model function as wrapped by the Model class for curve-fitting with an objective function for general purpose minimization with minimize or leastsq. An objective function will have a signature like:

     def objective(params, *args): 
    

    where params is a lmfit.Parameters instance that you would have to unpack within the function and *args are optional arguments. This function will return the array -- the objective -- that will be minimized in the least squares sense.

    In contrast, a model function used to create a curve-fitting Model will have a signature like

     def modelfunc(x, par1, par2, par3, ..., **kws):
    

    where x is the independent variable and par1 ... will contain the values for the model parameters. The model function will return an array that models the data.

    There are lots of examples using objective functions and using model functions at https://github.com/lmfit/lmfit-py/tree/master/examples. To model your data, you would do

    # define MODEL **not objective** function: returns the array to model the data
    def sinedecay(x, amp, shift, omega, decay):
        """Model a decaying sine wave""" . 
        return amp * np.sin(x*omega + shift) * np.exp(-x*x*decay)
    
    
    # create model:
    smodel = Model(sinedecay)
    
    # create a set of Parameters
    params = smodel.make_params(amp=10, decay=0.1, shift=0, omega=3)
    params.set('amp', min=0)
    params.set('shift',  min=-np.pi/2., max=np.pi/2)
    
    # do fit
    result = smodel.fit(data, params, x=x)