Problem fitting a piece-wise function with an exponential term

24 views
Skip to first unread message

spring21

unread,
Oct 20, 2018, 4:47:25 PM10/20/18
to lmfit-py
Hi,
    I am working with some data that I want to fit using a piecewise function that is one value before time = t' and at time t' it adds a constant term + an exponential decay.

So far I have been fitting using the following defined function:

def function_single(t, baseline, A, B, tau, tprime):
    result = baseline + np.heaviside(t-tprime, 1)*(A + B*np.exp(-(t-tprime)/tau))
    return result

and on most occasions it works fine, but occasionally I get an error:
ValueError: The input contains nan values

I think the problem might be that it runs into under/overflow issues, so I thought maybe it would be a good idea to make a composite model where I define a model like

Model = Exponential + StepModel

as these are the built-in functions, but I can't figure out how to make the exponential function have the form B*exp(|t-t'|/tau) where it is 0 for t-t' < 0.

Any suggestions or ideas on how to fit a function like this in a smarter way would be very welcome.
Thank you a lot!!

Matt Newville

unread,
Oct 21, 2018, 8:40:28 AM10/21/18
to lmfit-py
Proviving a minimial, complete example that shows the problem is always helpful.   Please provide one in the future.  Without seeing actual code, all I can do is speculate and give vague suggestions. 

If I had to guess, I would agree that the most likely culprit and source of the Nans or Infs is overflows from the exponentiation.   If the model function returns an array with Infs or Nans to the solver code, it is very hard to recover.  The most effective approach is to look for and recover from values that will generate Nans and Infs with your model function.  For your example model function, you may need to make sure that `tau` is not close to zero, and put an upper bound on the argument of the exponent, perhaps something like this:
 
def function_single(t, baseline, A, B, tau, tprime):
    arg = -(t-tprime) / max(1.e-200, tau)
    arg[np.where(arg> 200.)] = 200.0
    result = baseline + np.heaviside(t-tprime, 1)*(A + B*np.exp(arg)
    return result

That sort of check will avoid some sources of Infs.  Other common sources of nans are taking `sqrt()` and `log()` of negative numbers, and raising negative numbers to fractional powers.  

Hope that helps,

--Matt

Reply all
Reply to author
Forward
0 new messages