Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ModelResult.eval_uncertainty does not work for complex models #900

Closed
newville opened this issue Jul 3, 2023 · 2 comments
Closed

ModelResult.eval_uncertainty does not work for complex models #900

newville opened this issue Jul 3, 2023 · 2 comments

Comments

@newville
Copy link
Member

newville commented Jul 3, 2023

Description

ModelResult.eval_uncertainty is not working properly for models that generate complex data. Running the method gives a warning about casting complex to float, and returns a real array of the same length as the model.

See https://groups.google.com/g/lmfit-py/c/tqaLwSBH-Ts/m/HBDP0y4HEAAJ

A Minimal, Complete, and Verifiable example
import numpy as np
from lmfit import Model

def cmplx(f, omega, areal, aimag, off, sigma):
    return (areal*np.cos(f*omega + off) + 1j*aimag*np.sin(f*omega + off))*np.exp(-f/sigma)

f = np.linspace(0, 10, 501)
dat = cmplx(f, 4, 10, 5, 0.2, 4.5) + (0.1+0.2j)*np.random.normal(scale=0.25, size=len(f))


mod = Model(cmplx)
params= mod.make_params(omega=5, areal=5, aimag=5,
                        off={'value': 0.5, 'min': -2, 'max':2},
                        sigma={'value': 3, 'min': 1.e-5, 'max': 1000})

result = mod.fit(dat,  params=params, f=f)
print(result.fit_report())

dfit = result.eval_uncertainty()
print('Data  : ', len(dat), dat.dtype)
print('Fit   : ', len(result.best_fit), result.best_fit.dtype)
print('D_Fit : ', len(dfit), dfit.dtype)
Fit report:
[[Model]]
    Model(cmplx)
[[Fit Statistics]]
    # fitting method   = leastsq
    # function evals   = 49
    # data points      = 1002
    # variables        = 5
    chi-square         = 1.59641755
    reduced chi-square = 0.00160122
    Akaike info crit   = -6444.87518
    Bayesian info crit = -6420.32641
    R-squared          = (0.9996115675746108+2.8498175850731404e-05j)
[[Variables]]
    omega:  3.99961585 +/- 2.4279e-04 (0.01%) (init = 5)
    areal:  9.99808550 +/- 0.00755744 (0.08%) (init = 5)
    aimag:  4.99883013 +/- 0.00590546 (0.12%) (init = 5)
    off:    0.20069043 +/- 7.0226e-04 (0.35%) (init = 0.5)
    sigma:  4.50042395 +/- 0.00497052 (0.11%) (init = 3)
[[Correlations]] (unreported correlations are < 0.100)
    C(omega, off)   = -0.7328
    C(areal, sigma) = -0.6982
    C(aimag, sigma) = -0.4351
    C(areal, aimag) = +0.3021
/Users/Newville/Codes/lmfit-py/lmfit/model.py:1602: ComplexWarning: Casting complex values to real discards the imaginary part
  fjac[key][i] = (res1[key] - res2[key]) / (2*dval)
Data  :  501 complex128
Fit   :  501 complex128
D_Fit :  501 float64

Version information

Python: 3.10.10 | packaged by conda-forge | (main, Mar 24 2023, 20:17:34) [Clang 14.0.6 ]

lmfit: 1.2.1.post49+g78c9994, scipy: 1.11.0, numpy: 1.25.0,asteval: 0.9.30, uncertainties: 3.1.7

Link(s)

See https://groups.google.com/g/lmfit-py/c/tqaLwSBH-Ts/m/HBDP0y4HEAAJ

@newville
Copy link
Member Author

newville commented Jul 3, 2023

I have a proof-of-principle fix but will wait until #888 is merged to make a PR.

@newville
Copy link
Member Author

newville commented Jul 6, 2023

This is now fixed in the master branch, with a test added from the example code above.

@newville newville closed this as completed Jul 6, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant