Lmfit gives -1 correlation and large uncertainty (python)

Question

I am trying to fit a model function to a curve using the lmfit module.

The curve that I am fitting is set up as follows:

e(x) = exp(-(x-X)/x0) for x larger or equal than X, 0 otherwise.

G(x) = (1/sqrt(2*pi)*sigma) * exp(-x^2/2*sigma^2)

The model fit M(x) = E * conv(e,G)(x) + B

Where e is a truncated exponential, G is a gaussian, and E and B are constants. The operator between e and G is a convolution.

When I try to fit this function to my data I get a good fit. However, the fit is very sensitive to my initial value that I provide for X. This is also reflected in the uncertainty in the parameters:

[[Model]]
    ((Model(totemiss) * (Model(exptruncated) <function convolve at 0x7f139e2dcde8> Model(gaussian))) + Model(background))
[[Fit Statistics]]
    # fitting method   = leastsq
    # function evals   = 67
    # data points      = 54
    # variables        = 5
    chi-square         = 120558969110355112544642583094864038386991104.00000
    reduced chi-square = 2460387124701124853181382654239391973638144.00000
    Akaike info crit   = 5275.63336
    Bayesian info crit = 5285.57828
[[Variables]]
    E:           9.7316e+28 +/- 2.41e+33 (2475007.74%) (init= 1.2e+29)
    x0:          5.9420e+06 +/- 9.52e+04 (1.60%) (init= 5000000)
    X:           4.9049e+05 +/- 1.47e+11 (29978575.17%) (init= 100000)
    sigma:       2.6258e+06 +/- 5.74e+04 (2.19%) (init= 2000000)
    center:      0 (fixed)
    amplitude:   1 (fixed)
    B:           3.9017e+22 +/- 3.75e+20 (0.96%) (init= 4.5e+22)
[[Correlations]] (unreported correlations are <  0.100)
    C(E, X)                      = -1.000
    C(sigma, B)                  = -0.429
    C(x0, sigma)                 = -0.283
    C(x0, B)                     = -0.266
    C(E, x0)                     = -0.105
    C(x0, X)                     =  0.105

I suspect this has something to do due to the correlation between E and X being -1.00, which does not make any sense. I am trying to find out why I get this error and I believe it might be in the definition of the model:

    def exptruncated(x, x0, X):
        return np.exp(-(x-X)/x0)* (x > X)

        #Define convolution operator
    def convolve(arr, kernel):
        npts = min(len(arr), len(kernel))
        pad  = np.ones(npts)
        tmp  = np.concatenate((pad*arr[0], arr, pad*arr[-1]))
        out  = np.convolve(tmp, kernel, mode='valid')
        noff = int((len(out) - npts)/2)
        return out[noff:noff+npts]

       #Constant value for total emissions#
    def totemiss(x,E):
        return E

       #Constant value for background value
    def background(x,B):
        return B

        # create Composite Model using the custom convolution operator
        # M(x) = E + conv(exp,gauss) + B

    mod  = Model(totemiss)* CompositeModel(Model(exptruncated), Model(gaussian), convolve) + Model(background)
    mod.set_param_hint('x0',value=50*1e5,min=0,max=60*1e5)
    mod.set_param_hint('amplitude',value=1.0)
    mod.set_param_hint('center',value=0.0)
    mod.set_param_hint('sigma',value=20*1e5,min=0,max=100*1e5)
    mod.set_param_hint('X',value=1.0*1e5,min=0, max=5.0*1e5)
    mod.set_param_hint('B',value=0.45*1e23,min=0.3*1e23,max=1.0*1e23)
    mod.set_param_hint('E',value=1.2*1e29,min=1.2*1e26,max=1.0*1e32)


    pars = mod.make_params()

    pars['amplitude'].vary = False
    pars['center'].vary = False
    result =  mod.fit(y, params=pars, x=x)

    comps = result.eval_components(x=x)

Although I believe the model is the reason I am not able to find where the error comes from. Perhaps somebody of you can help me out!

[What a strange coincidence](https://stackoverflow.com/questions/49150218/fitting-curve-to-data-points-using-lmfit-python#comment85307778_49150218) Same course? — Mr. T, Mar 07 '18 at 22:01

score 0 · Answer 1 · answered Mar 07 '18 at 17:36

0

Why not just remove E from the model -- the X parameter is serving as a constant offset.

I'd also advise having parameters that are more reasonably scaled, closer to order of unity (roughly 1e-6 to 1e6, say). You can add scales of 1e10 and so on as needed in the model calculation, but it generally helps the calculations of covariance (used to determine how to update values in the fit) to have parameters more uniformly scaled.

answered Mar 07 '18 at 17:36

M Newville

7,486
2
16
29

X is an offset in the exponential, but the value of E (after the convolution is done) is just a scaling. How can it be that both of these are achieving the same effect? Isn't X just shifting the curve and E scaling the magnitude? – ksh Mar 08 '18 at 08:33
@ksh Well, did you plot the data? How many values do you have values with `x > X`? If `(x > X)` is doing nothing (I suspect this is the case), then `exp(X/x0)` is a constant scale factor to `exp(-x/x0)`, just as `E` is. `x0` interacts with the `x` array, but `X` probably does not. Using something like `x>X` rarely works for discrete data, but your values are so strangely scaled that I expect scaling is a larger problem. – M Newville Mar 08 '18 at 12:33

Lmfit gives -1 correlation and large uncertainty (python)

1 Answers1