snoplusuk / echidna

MIT License
4 stars 12 forks source link

Seperate floating and fixed backgrounds in limit setting #64

Closed jwaterfield closed 9 years ago

jwaterfield commented 9 years ago

Currently the code rescales all background spectra for each signal scaling. The code can be optimised by only scaling backgrounds which are being floated with the rest of the backgrounds scaled and fixed to their priors at the time of initialising the get_limit function.

I am currently thinking of separating backgrounds into 2 types which each get passed as a list. 1 where the backgrounds only have an accompanying prior and 1 where they have a prior and a numpy array of counts for floating the corresponding background. This will probably require quite a major rewrite of the limit_setting code. Do you agree before I precede @ashleyrback?

ashleyrback commented 9 years ago

Yes I think dividing them into two types of background, those that we want to scale and those that we don't seems an obvious solution, but I'm not sure if we actually need to pass _both_ lists to the SetLimit class. For example say we have something like:

# Load spectra
bkg1 = store.load("bkg1.hdf5")
bkg2 = store.load("bkg2.hdf5")

backgrounds = { bkg1: scaling1, bkg2: scaling2}  # realistically we probably want 
                                                 # this to be defined already in some
                                                 # sort of config file
# define backgrounds_to_float list of backgrounds that need to be floated

If you were passing this list to the SetLimit class, you'd then do:

# Create SetLimit instance
set_limit = limit_setting.LimitSetting(signal, backgrounds, backgrounds_to_float
                                       roi=roi, pre_shrink=True)

Then in SetLimit.__init__:

self._total_backgrounds = Spectra("total", total_events)
for background, scaling in backgrounds.iteritems():
    background.scale(scaling)
    self._total_backgrounds.add(background)

But you could equally do exactly the same thing before creating the SetLimit instance:

total_backgrounds = Spectra("total", total_events)
for background, scaling in backgrounds.iteritems():
    background.scale(scaling)
    total_backgrounds.add(background)

# Create SetLimit instance
set_limit = limit_setting.LimitSetting(signal, total_background, backgrounds_to_float
                                       roi=roi, pre_shrink=True)

here passing total_backgrounds as a pre-made background spectrum, for the backgrounds we don't want to float. Which I think avoids overcomplicating the constructor of SetLimit and since it is something that only needs to be done once. @jwaterfield let me know what you think.

Also I did have an alternate idea that we could take this opportunity to generalise the to any systematic (because when we float the backgrounds, what we are actually doing is floating the individual background rates as systematics). So what you could do is change the current _backgrounds list to a list of all _systematics and then create variations on the limit_config class for different types of systematic. But this would not be simple.

jwaterfield commented 9 years ago

Yeah that makes sense to me, it's similar to what I have coded currently although I have a function within LimitSetting which builds your fixed backgrounds for you rather than the user having to do this themselves but I guess it is more sensible to do this in a function outside the LimitSetting class.