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

Add detection utilities à la BgStats #569

Merged
merged 13 commits into from Jun 13, 2016
Merged

Add detection utilities à la BgStats #569

merged 13 commits into from Jun 13, 2016

Conversation

@jjlk
Copy link
Contributor

@jjlk jjlk commented Jun 9, 2016

Hi @cdeil @joleroi ,
I added a new class called ObservationStats in order to get statistics from an observation. The test functions are also here.

Julien

@cdeil cdeil added the feature label Jun 9, 2016
@cdeil cdeil added this to the 0.5 milestone Jun 9, 2016
@cdeil cdeil self-assigned this Jun 9, 2016
@cdeil cdeil changed the title Detection utilities à la BgStats Add detection utilities à la BgStats Jun 9, 2016
@@ -125,3 +126,33 @@ def show_in_browser(self):
"""Make HTML file and images in tmp dir, open in browser.
"""
raise NotImplementedError


# class ObservationSummary(object):

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

What's the plan with this class? Is it useful to keep this in a code comment here?
Maybe remove it and copy & paste it to somewhere else for yourself if you think it's useful?

This comment has been minimized.

@jjlk

jjlk Jun 12, 2016
Author Contributor

I try to implement this class and keep it simple. Now it take an ObservationStatsList and do some standard plots of detection. Maybe we could wrap it in an other class taking some options or a configuration file, to deal with the regions and the background substraction, as you suggested if I remembered correctly. We could also add the ObservationTableSummary utilities.

import astropy.units as u
from ...data import DataStore, ObservationList, ObservationStats, Target
from ...utils.testing import requires_data
from gammapy.extern.regions.shapes import CircleSkyRegion

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

Use relative imports for Gammapy here and on the following two lines.

@@ -63,7 +63,9 @@ def __str__(self):
keys = ['n_on', 'n_off', 'a_on', 'a_off',
'alpha', 'background', 'excess']
values = [self.n_on, self.n_off, self.a_on, self.a_off,
self.alpha, self.background(), self.excess()]
# self.alpha, self.background(), self.excess()]
### JLK

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

Thanks for fixing this. But please remove the ### JLK line.
We use git history to see who did what and never leave such comments in the code.

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

Since you've fixed up the Stats.__str__ method in this PR, would you mind adding a regression test that covers it as well?

return self.alpha_obs

@property
def sigma(self):

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

I think that property should be moved to the base class, no?
https://gammapy.readthedocs.io/en/latest/api/gammapy.stats.Stats.html

Note that it will still call the alpha from the subclass ... maybe you can just try moving it to the base class and see if all your tests still pass?

This comment has been minimized.

@jjlk

jjlk Jun 12, 2016
Author Contributor

In the class stats/Stats the property alpha returns the ratio of the exposures (a_on/a_off). When dealing with stacking of multiple observations, we usually average the normalisations with the number of off events, not with the mean acceptance. So I guess we need to re-implement the property alpha in the ObservationStats to deal with the weighting, right ?

This comment has been minimized.

@cdeil

cdeil Jun 13, 2016
Member

Ah yes, I forgot about that tricky part for the alpha computation.

Note that at some point I added this standalone function:
http://docs.gammapy.org/en/latest/api/gammapy.stats.combine_stats.html

This has some version of the alpha computation for two observations (probably not debugged / tested).
It should somehow be integrated with the Stats class better.
Maybe it's possible to have a stack classmethod on the base Stats class as well, i.e. follow the design pattern you implemented for ObservationStats?

So I guess we need to re-implement the property alpha in the ObservationStats to deal with the weighting, right ?

I'm not sure. Currently your design is such that the stacked ObservationStats doesn't contain arrays of n_on, a_on, ... but just the stacked values, right? Then the weighting has to happen at stacking time, not alpha computation time, no?

Is it clear to you how to handle / resolve this?
If no, I could have a closer look and think about it, or we delay this issue to a future pull request.
Let me know what you prefer ...

def _get_on_events(obs, target):
"""Number of ON events in the region of interest (`int`)
"""
print(target)

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

Remove print statement.

ss += 'Alpha: {}\n'.format(self.alpha)
ss += 'Bkg events in On region: {}\n'.format(self.background)
ss += 'Excess: {}\n'.format(self.excess)
ss += 'Gamma rate: {}\n'.format(self.n_on / self.livetime.to(u.min))

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

Maybe make the gamma and background rate properties on the class?
I think they are useful to have, and you're printing them, so it seems to be common info someone might want to get.

n_on : int
number of on events
n_off : int
number of on events

This comment has been minimized.

@cdeil

cdeil Jun 9, 2016
Member

There's a bunch of copy & paste errors ... all parameters are described as "number of on events"

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 9, 2016

I've left a bunch of minor comments on code cleanup.

I'll try to find some time tomorrow morning to think about the API and read through the implementation and tests in detail.

hijlk added 4 commits Jun 12, 2016
hijlk
…ots from an ObservationStatsList
@jjlk
Copy link
Contributor Author

@jjlk jjlk commented Jun 12, 2016

Hi @cdeil ,
I took into account your comments. I also implemented the class ObservationSummary and as I said above:

I tried to implement this class and keep it simple. Now it take an ObservationStatsList and do some standard plots of detection. Maybe we could wrap it in an other class taking some options or a configuration file, to deal with the regions and the background substraction, as you suggested if I remember correctly. We could also add the ObservationTableSummary utilities.

Maybe we could find some time do discuss this tomorrow via slack or Skype ?

Cheers !

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2016

There's two test fails on travis-ci:
https://travis-ci.org/gammapy/gammapy/jobs/137101665#L944
You should be able to reproduce those locally by running:

python setup.py test -V -t gammapy/data/tests/test_observation_stats.py

Fixing those is the most important part so that this pull request can go in.
Other things (stacked alpha computation, observation stats list to / from table, observation summary, ...) can be done in future pull requests one addition at a time.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2016

I see that you've implemented one version of ObservationStatsList to arrays in the ObservationSummary __init__ and _init_values methods.

I'm not sure if that is the best place for this code ... naively I would have expected methods ObservationStatsList.to_table and ObservationStatsList.from_table with an option on whether it's cumulative or not. But overall I'm not clear on where the distinction between ObservationStatsList and ObservationStats lies. I guess it's up to you @jjlk .
Eventually we also want a docs/data/stats.rst tutorial page that shows an example of how to use these classes.

But this PR is already 500 lines and hard for me to review, so I'm happy to merge as soon as tests pass and you think it's ready @jjlk , and then I'll use it for a bit (e.g. for the data release scripts) and then can give better feedback.

Thank you for this!

hijlk
@jjlk
Copy link
Contributor Author

@jjlk jjlk commented Jun 13, 2016

Hi @cdeil , I fixed the ObservationStats::str class and the test is fine now. I'll tackle your other comments in new pull requests

return ax

def __str__(self):
self.summary()

This comment has been minimized.

@cdeil

cdeil Jun 13, 2016
Member

There's one test that's still failing because you don't return the string here:
https://travis-ci.org/gammapy/gammapy/jobs/137196513#L944

Either put this:

def __str__(self):
    return self.summary()

or, if the summary method doesn't have any arguments to configure the summary, you can also get rid of it and just put it's implementation directly in __str__.

This comment has been minimized.

@jjlk

jjlk Jun 13, 2016
Author Contributor

@cdeil : Sorry about that, I didn't add the file. Now it's fixed

hijlk
"""Observation statistics.
Class allowing to summarize observation
(`~gammapy.data.DataStoreObservation) statistics

This comment has been minimized.

@cdeil

cdeil Jun 13, 2016
Member

One more fail on travis-ci in the docs build:
https://travis-ci.org/gammapy/gammapy/jobs/137211316#L1660

/home/travis/build/gammapy/gammapy/build/lib.linux-x86_64-2.7/gammapy/data/observation_stats.py:docstring of gammapy.data.ObservationStats:5: WARNING: Inline interpreted text or phrase reference start-string without end-string.

There's a missing backtick here, it should be:

(`~gammapy.data.DataStoreObservation`) statistics
run.pointing_radec, get_mask(), run.events)
stats = ObservationStats.from_target(run, get_target(), bg)
text = str(stats)
assert ('Observation summary report' in text)

This comment has been minimized.

@cdeil

cdeil Jun 13, 2016
Member

Can you please always remove these unused parentheses around text asserts?
(there's several places where you put this.)
I.e. instead of:

assert ('Observation summary report' in text)

write

assert 'Observation summary report' in text
return run


def get_target():

This comment has been minimized.

@cdeil

cdeil Jun 13, 2016
Member

Making this a pytest fixture would be a bit better.
https://github.com/jiffyclub/pytest-features#fixtures

I.e. put this here:

@pytest.fixture
def target():
    ... same as what you have now ...

and then inject target objects into the tests that use it, i.e. below:

def test_str(target):

At the moment we're not using fixtures a lot in Gammapy tests, but this is something that we'll gradually clean up ... almost always input test data will be injected into tests via fixtures.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2016

I promise I'll stop suggesting thing now on this PR and we'll merge as soon as travis-ci passes.
Thanks for your patience, @jjlk !

@jjlk
Copy link
Contributor Author

@jjlk jjlk commented Jun 13, 2016

Haha !
@cdeil : changes done

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2016

Looks good. I'm merging this in. Thanks!

@cdeil cdeil merged commit e1b31ad into gammapy:master Jun 13, 2016
0 of 2 checks passed
0 of 2 checks passed
continuous-integration/appveyor/pr Waiting for AppVeyor build to complete
Details
continuous-integration/travis-ci/pr The Travis CI build is in progress
Details
@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2016

@jjlk - I've done some minor cleanup in master (commits 25fe29d, 96e3cef, 110a741) for things I noticed when opening up the files locally in PyCharm. So please make sure to re-start from latest master if you continue editing these classes.

One more thing that PyCharm pointed out to me when I browsed your code locally is that the Parameters description for gammapy.data.ObservationStats.from_target is incorrect.
The method takes obs, target, and bg_estimate, and the parameter description says obs_table, target_pos and bg_method.
I didn't check what is actually the case or fix that one. @jjlk - Maybe you can fix this in your next PR. :-)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Linked issues

Successfully merging this pull request may close these issues.

None yet

2 participants
You can’t perform that action at this time.