Repository navigation
Gamma constraints 2 - #259
Conversation
for more information, see https://pre-commit.ci
for more information, see https://pre-commit.ci
Coverage Report for CI Build 34330474976Coverage decreased (-0.02%) to 74.708%Details
Uncovered Changes
Coverage Regressions1 previously-covered line in 1 file lost coverage.
Coverage Stats
💛 - Coveralls |
|
What's the status on this PR @Ananthu-Ravindran? 😊 |
|
Hey @hammannr , The PR is ready, I think. I haven't had time to do more with it since this note a few months ago. I wanted to do some checks on how the threshold function is affected by these constraints, but I am currently busy with the non-stats parts of the wimp analysis. |
|
I could probably have a look after my vacation. One thing I noticed is that unittests are still missing. Do you think you could add some? |
|
Sorry, reading this only now. I will add it to my todo list and try to do it this week. |
for more information, see https://pre-commit.ci
|
A week was about 2 months (sorry :-) ), but I think these tests should do? @hammannr Hmm, CodeFactor doesn't like my try except pass :/ |
Co-authored-by: Copilot <copilot@github.com>
Co-authored-by: Copilot <copilot@github.com>
Co-authored-by: Copilot <copilot@github.com>
Co-authored-by: Copilot <copilot@github.com>
hammannr
left a comment
There was a problem hiding this comment.
Thanks a lot @Ananthu-Ravindran for this nice PR! I now finally had the time to thoroughly check the code. I made a few minor changes and one change that affects how a change in n_sideband automatically changes the value of uncertainty. Please check if those changes make sense to you or otherwise feel free to revert them. I just thought it's more efficient to quickly do this myself.
There are still a few things that may need some attention to make sure everything works as expected. I think overall the PR is already in pretty good shape and I hope we can merge it soon! :)
| if not param.from_sideband: | ||
| parameter_meas = gen.rvs() | ||
| else: | ||
| parameter_meas = stats.poisson(mu=param.n_sideband * param.nominal_value).rvs() |
There was a problem hiding this comment.
Why do you multiply param.nominal_value here? Shouldn't this be simply the number of events in the sideband?
There was a problem hiding this comment.
At the moment this implementation doesn't take into account the generate values. So if you'd scan let's say the NR rate, 10 times higher NR background wouldn't reduce the uncertainty. I'm not sure if this is the intended behaviour. One could probably argue in both directions but it should be clearly stated in the docstring! No user should dig this info up in the code :D
There was a problem hiding this comment.
Or did you mean to use the central value (i.e. generate or nominal) one here? In any case, we should think of the consequences. Right now, if a user decides to use a nominal value of 10 they would get a narrower constraint than when they would absorb the factor 10 in the template and instead define the nominal value as 1. This can easily lead to trouble, I think.
There was a problem hiding this comment.
So after thinking a bit more, I think it will be almost impossible to avoid any ambiguity of the rescaling of template vs rescaling of the nominal value so I think the relative uncertainty should stay unaffected from this. So if some people use a nominal value that's not 1, the shape of the constraints shouldn't change, only their location should be rescaled to this nominal value. What do you think?
Of course, this then misses the fact that a larger background leads to better ancillary constraints. But this is also the case for all other uncertainties and the user would then need to manually do this by putting the change in rate and change in n_sideband in to_zip.
There was a problem hiding this comment.
Hey Robert, coming back to this now. Actually, the idea was that the uncertainty for constraints from sideband remain independent from the rate in the RoI. So, param.n_sideband is always the absolute number of events you expect to see in the sideband; and that alone dictates the width of the sideband.
But I think you are right, currently a larger rate multiplier will make the constraint width narrower. I am looking into this again; but I think it is not straightforward to do this. There may be a better way, but one way I can think of is by asserting that the param.nominal_value should always be 1 for sideband constraints. On the other hand, we can shift the location of the constraint using generate_values while maintaining its relative width solely based on param.n_sideband.
There was a problem hiding this comment.
Actually, I am not sure if this is also possible with discrete poisson toys. Let me check a little more
There was a problem hiding this comment.
I thinnk this is what I'd suggest: The param.n_sideband fully defines the width of the constraint. Everything else just scales the rate of the component but not the sideband i.e. the constraint width
There was a problem hiding this comment.
So, I think because we need the ancillary toys to be discrete and sampled from a Poisson, its not possible to fix the width of the constraint function when you scan the central values. The width of the poisson of count mu will have to be sqrt(mu).
If we set a fixed width and just move the same gamma constraint function up and down, then we set a hard wall for the constraint at mu < central_value-1. It also gives us situations where if we scan at toy_anc < 1, the constraint function goes into negative. This is not physical and I think not an option.
I think (I may be wrong and would be happy to be wrong in this case), if we need both a physical constraint function and discrete ancillary measurements, then the absolute width has to scale with the nominal/generate value where one is scanning. This will mean that when someone scans at higher value, they get a more narrower constraint (relative to the mean). The absolute width scales with sqrt(mu) and the relative by 1/sqrt(mu). This does again bring back the problem you raised about the ambiguity between template normalization and rate, but this has to be set at the config level. I think a n_sideband only makes sense when tied to the normalization of the template (such that the n_off/n_on ratio is always constant) and this change has to made at the config level.
I will now push a commit that will change the poisson sampling from Pois(n_sideband × nominal_value) to Pois(n_sideband × generate_value), so scanning either generate_value or nominal_value will give the exact same result in this type of constraint.
What do you think? If I'm missing something somewhere, please let me know :D
| # Expected sideband counts at rate multiplier = 1. Sets a Poisson (Gamma) constraint | ||
| # with relative width 1/sqrt(rate * n_sideband), so the width changes with the scanned | ||
| # rate; adjust n_sideband to keep it fixed or after rescaling the template normalization. |
There was a problem hiding this comment.
I'm still not convinced this is a behaviour users would expect. See e.g. the CES analysis, which uses ~number of events as the nominal rate values
No description provided.