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

Improve spectral model error propagation #2585

Merged
merged 21 commits into from Nov 21, 2019
Merged

Improve spectral model error propagation #2585

merged 21 commits into from Nov 21, 2019

Conversation

@cdeil
Copy link
Member

cdeil commented Nov 21, 2019

This PR changes the spectral model error propagation methods from the uncertainties-based implementation to a custom computation. The reasons for this change are described in detail in PIG 14 (see #2255).

The main change here is a new implementation for these methods:

  • SpectralModel.evaluate_error
  • SpectralModel.integral_error
  • SpectralModel.energy_flux_error
    A detailed test based on LogParabolaSpectralModel is added, to make sure the new implementation gives consistent results with the old implementation.

A regression test for the ECPL model case, closing #2007 now.

This PR should also fix #2190 (Naima models) and #1046 (absorbed models), since the new implementation is generic. I'm not adding regression tests for those here, but I've left a TODO comment that this should be done in the future.

These custom methods are deleted here:

  • PowerLawSpectralModel.integral_error
  • PowerLawSpectralModel.energy_flux_error
  • PowerLaw2SpectralModel.integral_error
  • PowerLaw2SpectralModel.energy_flux_error
    They are not useful, the implementation in the base class works for any model, including the power law variants. Performance is not a concern here, if at all, custom methods for evaluate_error should have been added in subclasses, because that's evaluated in ~ 100 points for spectral error band plots, wheres integral_error and energy_flux_error is usually evaluated for a single energy range.

Work in progres ...

@cdeil cdeil added this to the 0.15 milestone Nov 21, 2019
@cdeil cdeil requested a review from adonath Nov 21, 2019
@cdeil cdeil self-assigned this Nov 21, 2019
@cdeil cdeil added this to In progress in gammapy.modeling via automation Nov 21, 2019
@cdeil cdeil assigned adonath and unassigned cdeil Nov 21, 2019
@cdeil

This comment has been minimized.

Copy link
Member Author

cdeil commented Nov 21, 2019

@adonath - Could you please finish this PR? (or comment with suggestions)

There's two unit-related remaining fails that you can reproduce with:

pytest gammapy/catalog/tests/test_fermi.py

Probably the issue is the hard-coded energy.to_value("TeV") here?
https://github.com/gammapy/gammapy/pull/2585/files#diff-99eca274527206791c224179cd06b3aeR59

Another possible unit issue is this _convert_energy here, which relies on hard-coded parameters names (which is error prone, and will break soon when we change if not already broken for that Fermi case):
https://github.com/gammapy/gammapy/pull/2585/files#diff-99eca274527206791c224179cd06b3aeR46-R52

Finally, this should be improved:

# This function is copied over from https://github.com/zblz/naima/blob/master/naima/utils.py#L261
# and slightly modified to allow use with the uncertainties package

This can be deleted for sure:

# arrays with uncertainties contain objects
if y.dtype == "O":
from uncertainties.unumpy import log10
# uncertainties.unumpy.log10 can't deal with tiny values see
# https://github.com/gammapy/gammapy/issues/687, so we filter out the values
# here. As the values are so small it doesn't affect the final result.
# the sqrt is taken to create a margin, because of the later division
# y[slice2] / y[slice1]
valid = y > np.sqrt(np.finfo(float).tiny)
x, y = x[valid], y[valid]

but probably there's other improvements that can be made in that function. E.g. below the same slice is taken multiple times which is slow. If we use integrate this will be in our model fitting hot loop, no?
Do we have a strategy for unit handling in the helper integrate and gradient methods?
Maybe it's easiest to keep units out of those completely?

Finally, I removed integral_error and energy_flux_error for now.
Do we need / want it back?
If yes, it would have to be added back using the new implementation.
I'm OK to just remove them, but if it's easy to bring back, OK.
There's already an xfailing test that could be used for those in test_spectrum.py.

@adonath

This comment has been minimized.

Copy link
Member

adonath commented Nov 21, 2019

@cdeil I fixed the remaining test fails and did a bit of clean-up.

adonath and others added 2 commits Nov 21, 2019
Copy link
Member

adonath left a comment

Thanks @cdeil! I have no further comments...

@adonath adonath merged commit cdaedda into gammapy:master Nov 21, 2019
7 of 10 checks passed
7 of 10 checks passed
greeting
Details
Scrutinizer Running
Details
continuous-integration/travis-ci/pr The Travis CI build is in progress
Details
gammapy.gammapy in progress
Details
Codacy/PR Quality Review Up to standards. A positive pull request.
Details
gammapy.gammapy (DevDocs) DevDocs succeeded
Details
gammapy.gammapy (Lint) Lint succeeded
Details
gammapy.gammapy (Test Python36) Test Python36 succeeded
Details
gammapy.gammapy (Test Windows36) Test Windows36 succeeded
Details
gammapy.gammapy (Test Windows37) Test Windows37 succeeded
Details
gammapy.modeling automation moved this from In progress to Done Nov 21, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
2 participants
You can’t perform that action at this time.