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

Uncertainties for constrained parameters #1682

Merged
merged 10 commits into from
Aug 31, 2021
Merged

Uncertainties for constrained parameters #1682

merged 10 commits into from
Aug 31, 2021

Conversation

m2cci-NMZ
Copy link
Contributor

This PR aims at having uncertainties computed for constrained parameters and displayed in the GUI. The errors are popagated using the uncertainties package, as suggested by @pkienzle during the last biweekly call.

m2cci-NMZ and others added 3 commits September 30, 2020 13:07
Update fitted parameters values with uncertainty objects and propagate them to constrained parameters
Add uncertainty objects with 0 standard deviation to non fitted parameters. This makes sure that parameter value attribute is always an uncertainty object
Copy link
Contributor

@pkienzle pkienzle left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Code looks fine. We should probably pin a specific version of the package before testing and release.

@RichardHeenan
Copy link
Contributor

Presume someone has to add the dependency package here so it gets into release & developer code.

@RichardHeenan
Copy link
Contributor

RichardHeenan commented Oct 6, 2020

Have now done pip install uncertainties, so can build locally on Windows.
When I run my local build I cannot successfully enter any constraints, so have not yet managed to test the new feature.. SORTED now - had not assigned model to thrid FitPage, so tripped an older bug that is fixed in the constraint branch. Might be good to merge this branch with the constraints one for a final test? Will do some testing now ....

See log below.

11:04:06 - ERROR: Traceback (most recent call last):
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\ConstraintWidget.py", line 945, in onAcceptConstraint
constrained_tab.addConstraintToRow(constraint, constrained_row)
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FittingWidget.py", line 810, in addConstraintToRow
symbol_dict = self.parent.perspective().getSymbolDictForConstraints()
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FittingPerspective.py", line 537, in getSymbolDictForConstraints
symbol_dict.update(tab.getSymbolDict())
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FittingWidget.py", line 4375, in getSymbolDict
model_name = self.kernel_module.name
AttributeError: 'NoneType' object has no attribute 'name'

@RichardHeenan
Copy link
Contributor

RichardHeenan commented Oct 6, 2020

So now I am slightly puzzled, but may not be clever enough, in my local build I did
M1:sld_shell = 0.6 *M2.sld_core+0.4*M2.sld_solvent (sorry it changed asterix to italic)
In M2, sld_core = 6.5025 (0.0094271)
sld_solvent= 10.893 (0.0072262)

But M1.sld_shell = 8.2588 (0.006352)

I was expecting a simple linear combination of the sigma here,
8.2947 (0.6*0.0094 + 0.4*0.00722) = 8.2947 (0.0085)
and the actual value of the param is wrong, so is something out of step inthe updating or passing of parameter values ?
On another, slightly different, fit I did get the correct value , but stil the sigma is smaller than that of the two parents. Is is allowing for parameter correlation?

@m2cci-NMZ
Copy link
Contributor Author

So now I am slightly puzzled, but may not be clever enough, in my local build I did
M1:sld_shell = 0.6 M2.sld_core+0.4M2.sld_solvent (sorry it changed asterix to italic)
In M2, sld_core = 6.5025 (0.0094271)
sld_solvent= 10.893 (0.0072262)

But M1.sld_shell = 8.2588 (0.006352)

I was expecting a simple linear combination of the sigma here,
8.2947 (0.60.0094 + 0.40.00722) = 8.2947 (0.0085)
and the actual value of the param is wrong, so is something out of step inthe updating or passing of parameter values ?

It is a linear combination for variance, so it should be sqrt((0.6*0.0094)**2 + (0.4*0.00722)**2) :-)

On another, slightly different, fit I did get the correct value , but stil the sigma is smaller than that of the two parents. Is is allowing for parameter correlation?

Parameter correlation is not taken into account for the uncertainty calculation. It thus assumes parameters are independent.

@RichardHeenan
Copy link
Contributor

In which case I am concerned that the package is not doing what most people expect for independent variables.

(I did check, that as its docs say M1.sld_solvent = M2.sld_solvent - M2.sld_solvent gives zero with error zero. Which is also possibly not good., unless it is actaully tracking the algebra properly. )

@pkienzle
Copy link
Contributor

pkienzle commented Oct 6, 2020

Re: formating, use backticks around code and triple-backticks around code blocks. This will keep asterisks as *.

Re: results:

>>> 0.6*6.5025 + 0.4*10.893                                                                                                                                                                                                                                                
8.258700000000001

I don't know how you are calculating 8.2947.

Uncertainty via forward Monte Carlo matches the results of the uncertainty package:

>>> core, solv = U(6.5025, 0.0094271), U(10.893, 0.00722) 
>>> 0.6*core + 0.4*solv                                                                                                                                                                                                                
8.258700000000001+/-0.006350891369532312
>>> A, B = np.random.normal(6.5025, 0.0094271, size=1000000), np.random.normal(10.893, 0.00722, size=1000000)

>>> std(A), std(B)
(0.009417652426904453, 0.007223868793688441)

>>> mean(0.6*A+0.4*B), std(0.6*A+0.4*B, ddof=1)
(8.258705161496879, 0.006342674075400538)

It differs from SasView in the fourth digit , presumably because of limited precision in the intermediate printed representation.

This also matches the analytical variance as given on wikipedia:

>>> sqrt((0.6*0.0094271)**2 + (0.4*0.00722)**2)
0.006350891369532312

It is actually tracking the algebra properly, so (A+B)/A will take into account the covariance between numerator and denominator in the uncertainty propagation:

>>> mean((A+B)/A), std((A+B)/A, ddof=1)
(2.6752054928106235, 0.0026700115120890005)

>>> (core + solv)/core
2.67520184544406+/-0.002670430600778083

@pkienzle
Copy link
Contributor

pkienzle commented Oct 6, 2020

As Nicholas mentioned, the code does not encode the parameter covariance.

It turns out that this would be pretty easy to do: https://pythonhosted.org/uncertainties/user_guide.html#use-of-a-covariance-matrix

With dream, the covariance could be computed directly from the sample instead of using the Hessian approximation.

@RichardHeenan
Copy link
Contributor

OK, they didn't teach it like that back in my "high school" days. I suppose the simple view is that combining the two observations results in a smaller fractional uncertainty as there are now more observations. I am happy now, so will approve.

If we could put in the covariances also, that would be very good.

[I don't know how you are calculating 8.2947 - sorry must have mistyped in my calculator, at least twice over]

Take into account parameter correlation when computing the uncertainties
@m2cci-NMZ
Copy link
Contributor Author

I've updated the code so parameter covariance is taken into account when computing the uncertainties

# TODO: move uncertainty propagation into bumps
# TODO: should scale stderr by sqrt(chisq/DOF) if dy is unknown
values, errs, state = result['value'], result["stderr"], result[
'uncertainty']
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The call to problem.cov must have the returned fit values inserted into the model. Need to guarantee that is the case even if the fit is being run in a different thread. It will be safer to add cov to the result object at the same time that stderr was added since we know the model is in the correct state.

Then we can write:

from uncertainties import ufloat, correlated_values
from bumps.parameters import unique as unique_parameters
...
# Turn all parameters (fixed and varying) into uncertainty objects with zero uncertainty.
for p in unique_parameters(problem.model_parameters()):
    p.value = ufloat(p.value, 0)
# Update varying parameters with fitted uncertainty, preserving correlations.
fitted = correlated_values(result['value'], result['cov'])
for p, v in zip(self._parameters, fitted):
    p.value = v
# Propagate correlated uncertainty through constraints.
problem.setp_hook()

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I still have to separate the case where only one parameter is being fitted since np.cov returns a 0 dimension array when a (1,n) array is passed as argument. This will cause correlated_values to fail since it is excepting at least a 1-d array.

Have the result dict return the covariance matrix to make sure that the probem is in the correct state when computing the matrix
@RichardHeenan
Copy link
Contributor

In this branch today, fitting 3 sets simultaneously, but NO constraints anywhere. I am seeing errors as below, will try to work out what kicked them off...

12:16:39 - INFO: 2020-10-13 12:16:39
12:16:44 - INFO: 2020-10-13 12:16:44
12:28:12 - INFO: 2020-10-13 12:28:12
12:28:12 - ERROR: Fitting failed: Traceback (most recent call last):
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 79, in compute
result = list(map(map_apply, inputs))
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 19, in map_apply
return arguments0
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 16, in map_getattr
return getattr(classInstance, classFunc)(*args)
File "C:\sasview42\sasview\src\sas\sascalc\fit\BumpsFitting.py", line 308, in fit
fitted = correlated_values(values, cov)
File "C:\anaconda\envs\qt5\lib\site-packages\uncertainties\core.py", line 186, in correlated_values
tags)
File "C:\anaconda\envs\qt5\lib\site-packages\uncertainties\core.py", line 230, in correlated_values_norm
(variances, transform) = numpy.linalg.eigh(correlation_mat)
File "C:\anaconda\envs\qt5\lib\site-packages\numpy\linalg\linalg.py", line 1410, in eigh
w, vt = gufunc(a, signature=signature, extobj=extobj)
File "C:\anaconda\envs\qt5\lib\site-packages\numpy\linalg\linalg.py", line 95, in _raise_linalgerror_eigenvalues_nonconvergence
raise LinAlgError("Eigenvalues did not converge")
numpy.linalg.linalg.LinAlgError: Eigenvalues did not converge

12:28:12 - ERROR: Eigenvalues did not converge
None

12:28:12 - ERROR: Traceback (most recent call last):
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FittingWidget.py", line 1796, in fitComplete
self.chi2 = res.fitness
AttributeError: 'str' object has no attribute 'fitness'

@RichardHeenan
Copy link
Contributor

RichardHeenan commented Oct 13, 2020

Loaded & set up 3 fitpages, did a simultaneous fit OK, then a fit on fitpage1 OK, then another simultaneous fit - Error, slightly different to the one above:
12:39:54 - INFO: 2020-10-13 12:39:54
12:39:54 - ERROR: Fitting failed: Traceback (most recent call last):
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 79, in compute
result = list(map(map_apply, inputs))
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 19, in map_apply
return arguments0
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 16, in map_getattr
return getattr(classInstance, classFunc)(*args)
File "C:\sasview42\sasview\src\sas\sascalc\fit\BumpsFitting.py", line 308, in fit
fitted = correlated_values(values, cov)
File "C:\anaconda\envs\qt5\lib\site-packages\uncertainties\core.py", line 186, in correlated_values
tags)
File "C:\anaconda\envs\qt5\lib\site-packages\uncertainties\core.py", line 230, in correlated_values_norm
(variances, transform) = numpy.linalg.eigh(correlation_mat)
File "C:\anaconda\envs\qt5\lib\site-packages\numpy\linalg\linalg.py", line 1410, in eigh
w, vt = gufunc(a, signature=signature, extobj=extobj)
File "C:\anaconda\envs\qt5\lib\site-packages\numpy\linalg\linalg.py", line 95, in _raise_linalgerror_eigenvalues_nonconvergence
raise LinAlgError("Eigenvalues did not converge")
numpy.linalg.linalg.LinAlgError: Eigenvalues did not converge

12:39:54 - ERROR: Eigenvalues did not converge
None

Some further single fits on the three fitpages then work, some give errors. (Also why is residuals plot jumping up from where I parked it out of the way.)

@m2cci-NMZ
Copy link
Contributor Author

Loaded & set up 3 fitpages, did a simultaneous fit OK, then a fit on fitpage1 OK, then another simultaneous fit - Error, slightly different to the one above:
12:39:54 - INFO: 2020-10-13 12:39:54
12:39:54 - ERROR: Fitting failed: Traceback (most recent call last):
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 79, in compute
result = list(map(map_apply, inputs))
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 19, in map_apply
return arguments0
File "C:\sasview42\sasview\src\sas\qtgui\Perspectives\Fitting\FitThread.py", line 16, in map_getattr
return getattr(classInstance, classFunc)(*args)
File "C:\sasview42\sasview\src\sas\sascalc\fit\BumpsFitting.py", line 308, in fit
fitted = correlated_values(values, cov)
File "C:\anaconda\envs\qt5\lib\site-packages\uncertainties\core.py", line 186, in correlated_values
tags)
File "C:\anaconda\envs\qt5\lib\site-packages\uncertainties\core.py", line 230, in correlated_values_norm
(variances, transform) = numpy.linalg.eigh(correlation_mat)
File "C:\anaconda\envs\qt5\lib\site-packages\numpy\linalg\linalg.py", line 1410, in eigh
w, vt = gufunc(a, signature=signature, extobj=extobj)
File "C:\anaconda\envs\qt5\lib\site-packages\numpy\linalg\linalg.py", line 95, in _raise_linalgerror_eigenvalues_nonconvergence
raise LinAlgError("Eigenvalues did not converge")
numpy.linalg.linalg.LinAlgError: Eigenvalues did not converge

12:39:54 - ERROR: Eigenvalues did not converge
None

Some further single fits on the three fitpages then work, some give errors. (Also why is residuals plot jumping up from where I parked it out of the way.)

Could you send me the project causing the error? I can't figure out if the problem originates from some malformed covariance matrix or if it's something else...

@pkienzle
Copy link
Contributor

pkienzle commented Oct 14, 2020

It may depend on the particular version of LAPACK installed:

That'll make it kind of hard to debug. I tried to reproduce it with an ill-conditioned matrix and with a non-symmetric matrix, but no luck. Even a matrix of ones didn't throw an error.

If you print out the covariance matrix on error by changing this code to (untested):

import logging
try:
    fitted = correlated_values(values, cov)
except:
    fitted = [ufloat(v, dv) for v, dv in zip(values, errs)]
    with np.printoptions(precision=15): 
        logging.error(f"Parameter correlations not propagated. Covariance failed with cov=\n{cov}")

then we can at least see what kind of matrices cause problems.

If we only do this if there are derived parameters, then we won't flood the sasview log file, but it'll take longer to debug.

@m2cci-NMZ
Copy link
Contributor Author

Problem seems to come from negative values in the covariance diagonal which then fill the correlation matrix corresponding column and row with nan. I could sourround the failing function with try...except blocks and compute uncertainties without covariance in case it fails and warn the UI so the user gets some feedback.

@RichardHeenan
Copy link
Contributor

In my cases that gave the errors there were not any constraints set up, I just used "fit" on the simultaneous fit tab to run three fitpages at once. So perhaps code should not be trying to use the covariances in such circumstance?

@wpotrzebowski wpotrzebowski changed the base branch from ESS_GUI to main October 26, 2020 07:21
@butlerpd
Copy link
Member

There are a couple of "solutions" but need to do something as it is used in all fits, not just constrained.

@butlerpd
Copy link
Member

@pkienzle suggests using diagonal rather than NaN as an option. Agreed.

@RichardHeenan
Copy link
Contributor

RichardHeenan commented Oct 27, 2020

In this one try "fit" on fitpage2
cov_errors_failing2json.txt

this one has too many params fitting, so there will be extreme correlation. The fit engine does not always come up with the error, so depends what route the fit takes to get there.

In this one, "fit" on Simult Fit tab
cov_errors_failingjson.txt
here fitpage1 is the issue, thickness has gone small, but then setting it to say 12, and hit fit on that page, sasview hangs in the fit.

In my local build on Win10 of ESS_GUI_uncertainties both will give errors.

(I am presuming that the L-M engine is clever enough to re-invert the least squares matrix, without the Marquardt modification to the diagonal, before extracting the parameter variances etc, else if the fit has "converged" with a large marquardt lambda value, poor results (too small?) will appear for the variances & covariances.)

@RichardHeenan
Copy link
Contributor

Some related issues:
#628 singular covariance matrices
#994 show users the correlation matrix, to which I would add that even if this is tricky, can we warn users that the fit is close to singular or ill conditioned (too many parameters, highly correlated parameters, or parameters that have no effect at all on the fit)
#571 add values of correlation coeffs to DREAM plots.

Check if constraints are applied before calling `correllated_values`
Set covariance to `np.nan` if covariance computation fails
Ignore covariance in case `correlated_values` call raises `LinALgError` and disable error propagation in case it raises a `TypeError` (which is raised when using `math` functions).
@Caddy-Jones
Copy link
Contributor

Caddy-Jones commented Aug 10, 2021

Hey guys I've started looking at this PR.
I'm attempting to test the code, but cannot load the cov_errors_failing2json.txt or the cov_errors_failingjson.txt file into SasView.

This is the error message in receiving:
"When contacting the SasView team, mention the following:
Error: None
while loading Data:
cov_errors_failing2json.txt"

I've uploaded the files in txt format, and tried converting them to json.

@smk78
Copy link
Contributor

smk78 commented Aug 10, 2021

Hi @Caddy-Jones , have you tried renaming the file to .json (the .txt may have been added to add it to Github) and loading it as an analysis or project from the File menu? Just a thought...

@Caddy-Jones
Copy link
Contributor

That works, thanks :)

@Caddy-Jones
Copy link
Contributor

Hi @RichardHeenan, are you still experiencing the same errors with cov_errors_failingjson.txt and cov_errors_failing2json.txt ? Because I am not getting any error messages when running the fits for either of these files. I am wondering if the commits by the previous intern have fixed the problem.

@butlerpd butlerpd marked this pull request as ready for review August 17, 2021 13:26
@RichardHeenan
Copy link
Contributor

Hi @RichardHeenan, are you still experiencing the same errors with cov_errors_failingjson.txt and cov_errors_failing2json.txt ? Because I am not getting any error messages when running the fits for either of these files. I am wondering if the commits by the previous intern have fixed the problem.

hi @Caddy-Jones , when I load test projects as above into a local build of this branch, I no longer get error messages, and some uncertainties are appearing on the constrained parameters, so that is good!

The project loader is alas not overwriting Fitpage1, but adds FitPage2 through FitPage4

(Later I will try some freshly generated constrained fits and see what I can break, meanwhile you may as well merge all this and then we see if anyone notices the changes)

@butlerpd
Copy link
Member

If I understand correctly then @RichardHeenan , we are happy to merge this now but there may be some other issues that need to be ticketed separately? It looks like @Caddy-Jones has fixed the conflicts already so it is ready to merge otherwise.

@krzywon krzywon merged commit c5467c0 into main Aug 31, 2021
Caddy-Jones added a commit that referenced this pull request Sep 17, 2021
butlerpd added a commit that referenced this pull request Oct 12, 2021
@krzywon krzywon deleted the ESS_GUI_uncertainties branch August 18, 2022 14:52
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

Successfully merging this pull request may close these issues.

7 participants