[SciPy-User] curve_fit - fitting a sum of 'functions'

[Meesters, Christian]

Hi,

I'm dealing with a sigmoid curve which is pretty noisy at its maximum and the maximum is about three orders of magnitude bigger than the minimum. Just by eye I see that the data contain more than one 'mode' or, in other words, these curves have more than one step. In order to describe the data I would like to fit a sum normal cdf's like:

def normald(xdata, *args):
    if len(args) == 1: args = args[0]
    if len(args) < 4 or (len(args) - 1) % 3 != 0:
        raise AssertionError("Argument list too short or incomplete")
    rec        = np.zeros(len(xdata))
    base      = args[0]
    for amp, e50, var in zip(args[1::3], args[2::3], args[3::3]):
        rec += norm.cdf(xdata, loc = e50, scale = var) * amp
    return rec + base

(If somebody knows a more elegant way to describe this, I would be glad to read it. Fitting common exponential works as good or bad and the data is readily interpreted with cdfs of the normal distribution.)

Anyway, this than is passed to scipy's curve_fit:

from scipy.optimize import curve_fit

popt, pcov = curve_fit(normald, xdata, ydata, p0 = guess)

My problem now is with the 'guess' parameter: If handing over parameters as estimated by eye, all works fine, but more general approaches failed so far. I tried several approaches where the 'guesses' for the additional modes were appended within a loop and subsequently fitted, e.g. 'substract the fit from the original data and guess maximum, minimum, inflection point and append the new parameters to 'guess''. This, of course, is naive, because the upper part is so noisy and because the maximum is so much higher than the other modes, an least squares fit first tries to fit the noise. Truncating the data above the highest inflection point has no effect, so: Does anybody have an idea how to approach this problem in a more elegant / robust way?

Any pointer is appreciated. (Do I need to describe the problem better?)

TIA,
Christian

[Yury V. Zaytsev]

On Wed, 2011-05-25 at 07:36 +0000, Meesters, Christian wrote:
> My problem now is with the 'guess' parameter: If handing over
> parameters as estimated by eye, all works fine, but more general
> approaches failed so far.

I don't know if it would seem to make any sense to you, but you could
heavily filter your function with something like Savitzky–Golay
smoothing filter (it's basically just a weighted sum, so no numerical
heavy-lifting is involved) and use the filtered version to derive your
guesses, if the noise is really the problem...

--
Sincerely yours,
Yury V. Zaytsev

_______________________________________________
SciPy-User mailing list
SciPy...@scipy.org
http://mail.scipy.org/mailman/listinfo/scipy-user

[Meesters, Christian]

Oh, thanks. I forgot to mention that I tried a) a Savitzky-Golay filter, b) Wiener-Filtering and c) B-Spline smoothing. No effect, then I stopped. Point is that the noise is just a problem on the plateau and both, filtering and smoothing, either level this noise and part of the actual data or - depending on the parameters - follow and enhance the noise.

Yet, the noise can be truncated and still the fits aren't stable.
________________________________________
From: scipy-use...@scipy.org [scipy-use...@scipy.org] on behalf of Yury V. Zaytsev [yu...@shurup.com]
Sent: Wednesday, May 25, 2011 9:43 AM
To: SciPy Users List
Subject: Re: [SciPy-User] curve_fit - fitting a sum of 'functions'

[David Baddeley]

Have you got any way of determining the magnitude of the noise? You could try a
weighted fit (although you might have to take a step backwards and use lsqnonlin
- I'm not sure if curve_fit deals with weights. I think the spline fitting
functions can also take weights, which might make the spline smoothing approach
worth another shot.