Benjamin
Benjamin

Reputation: 11860

Dynamically define a function

I am trying to write a curve fitting function which returns the optimal parameters a, b and c, here is a simplified example:

import numpy
import scipy
from scipy.optimize import curve_fit

def f(x, a, b, c):
    return x * 2*a + 4*b - 5*c

xdata = numpy.array([1,3,6,8,10])
ydata = numpy.array([  0.91589774,   4.91589774,  10.91589774,  14.91589774,  18.91589774])
popt, pcov = scipy.optimize.curve_fit(f, xdata, ydata)

This works fine, but I want to give the user a chance to supply some (or none) of the parameters a, b or c, in which case they should be treated as constants and not estimated. How can I write f so that it fits only the parameters not supplied by the user?

Basically, I need to define f dynamically with the correct arguments. For instance if a was known by the user, f becomes:

def f(x, b, c):
    a = global_version_of_a
    return x * 2*a + 4*b - 5*c

Upvotes: 8

Views: 2120

Answers (5)

deprecated
deprecated

Reputation: 2185

There is already a package that does this:

https://lmfit.github.io/lmfit-py/index.html

From the README:

"LMfit-py provides a Least-Squares Minimization routine and class with a simple, flexible approach to parameterizing a model for fitting to data. Named Parameters can be held fixed or freely adjusted in the fit, or held between lower and upper bounds. In addition, parameters can be constrained as a simple mathematical expression of other Parameters."

Upvotes: 0

robince
robince

Reputation: 10967

I usually use a lambda for this purpose.

user_b, user_c = get_user_vals()
opt_fun = lambda x, a: f(x, a, user_b, user_c)
popt, pcov = scipy.optimize.curve_fit(opt_fun, xdata, ydata)

Upvotes: 1

unutbu
unutbu

Reputation: 879391

Taking a page from the collections.namedtuple playbook, you can use exec to "dynamically" define func:

import numpy as np
import scipy.optimize as optimize
import textwrap

funcstr=textwrap.dedent('''\
def func(x, {p}):
    return x * 2*a + 4*b - 5*c
''')
def make_model(**kwargs):
    params=set(('a','b','c')).difference(kwargs.keys())
    exec funcstr.format(p=','.join(params)) in kwargs
    return kwargs['func']

func=make_model(a=3, b=1)

xdata = np.array([1,3,6,8,10])
ydata = np.array([  0.91589774,   4.91589774,  10.91589774,  14.91589774,  18.91589774])
popt, pcov = optimize.curve_fit(func, xdata, ydata)
print(popt)
# [ 5.49682045]

Note the line

func=make_model(a=3, b=1)

You can pass whatever parameters you like to make_model. The parameters you pass to make_model become fixed constants in func. Whatever parameters remain become free parameters that optimize.curve_fit will try to fit.

For example, above, a=3 and b=1 become fixed constants in func. Actually, the exec statement places them in func's global namespace. func is thus defined as a function of x and the single parameter c. Note the return value for popt is an array of length 1 corresponding to the remaining free parameter c.


Regarding textwrap.dedent: In the above example, the call to textwrap.dedent is unnecessary. But in a "real-life" script, where funcstr is defined inside a function or at a deeper indentation level, textwrap.dedent allows you to write

def foo():
    funcstr=textwrap.dedent('''\
        def func(x, {p}):
            return x * 2*a + 4*b - 5*c
        ''')

instead of the visually unappealing

def foo():
    funcstr='''\
def func(x, {p}):
    return x * 2*a + 4*b - 5*c
'''

Some people prefer

def foo():
    funcstr=(
        'def func(x, {p}):\n'
        '    return x * 2*a + 4*b - 5*c'
        )

but I find quoting each line separately and adding explicit EOL characters a bit onerous. It does save you a function call however.

Upvotes: 8

Jan
Jan

Reputation: 735

If you want a simple solution based on curve_fit, I'd suggest that you wrap your function in a class. Minimal example:

import numpy
from scipy.optimize import curve_fit

class FitModel(object):
    def f(self, x, a, b, c):
        return x * 2*a + 4*b - 5*c

    def f_a(self, x, b, c):
        return self.f(x, self.a, b, c)

# user supplies a = 1.0
fitModel = FitModel()
fitModel.a = 1.0

xdata = numpy.array([1,3,6,8,10])
ydata = numpy.array([  0.91589774,   4.91589774,  10.91589774,  14.91589774,  18.91589774])
initial = (1.0,2.0)
popt, pconv = curve_fit(fitModel.f_a, xdata, ydata, initial)

Upvotes: 0

urschrei
urschrei

Reputation: 26859

def f(x, a = 10, b = 15, c = 25):
    return x * 2*a + 4*b - 5*c

If the user doesn't supply an argument for the parameter in question, whatever you specified on the right-hand side of the = sign will be used:

e.g.:
f(5, b = 20) will evaluate to return 5 * 2*10 + 4*20 - 5*25 and
f(7) will evaluate to return 7 * 2*10 + 4*15 - 5*25

Upvotes: -2

Related Questions