Non continuous limits

[jmsardain]

Hello,
We are using pyhf v0.7.0, and we are trying to set some limits. However we are getting a non continuous behavior of the observed CLs. We are using the function pyhf.infer.intervals.upperlimit to get the limits.

We have tried to vary different parameters of the optimizer (using pyhf.optimize.minuit_optimizer). We tried tolerance values of {0.02, 0.001, 0.0001}, strategy=0,1,2, steps=1000. However, these jumps still appear.
We also tried adding an offset to the negative log-likelihood bringing it closer to O(0-1) to make the optimisation more robust (we suspected numerical precision issues causing the non-continuous observed CLs). However even with this adjustment, we still see the problem.

Any guidance on how to make this converge would be appreciated.

Tagging @lhenkelm

[matthewfeickert]

@jmsardain I saw that this was also posted to the ATLAS Statistics Forum mailing list, so there will be some cross-talk, but I'm glad that you also posted this here so that there can be a public discussion even though the workspace itself won't be.

As @lukasheinrich has already mentioned on the mailing list

This kink is an expected feature and relates to the fact that there is an ambiguity how to reconstruct the best fit value from a 0 test statistic (since all best for values above the test value get mapped to zero for the exclusion test statistics).

Sorry for being brief, I can explain in detail once I’m on my laptop but maybe it already helps to know that this is expected.

As @will-cern mentioned, if you privately share your workspace (e.g. through a private GitHub Gist) we can walk through this a bit, though I think @lukasheinrich will also follow up here with a pedagogical example.

(perhaps this would be a good example for us to include in the pyhf users guide Jupyter Book, though there is already a start on this as a "learn" notebook from @lukasheinrich in PR #1476, but is requires more explanation to be instructionally useful).

[kratsg]

This is discussed a bit in #993 as well. You can add calc_kwargs = {'base_dir': 'clipped_normal'} to the calculator in order to get clipping as described in the PR.

[lhenkelm]

Thank you for pointing out the test statistic differences! We did use the default/non-clipped variant in the plot we shared above,
and we since switched to using the clipped version for our CLs scans.

However there is an issue around "jumpy" or non-monotonous p-values apart from the question of the test statistic you covered.
I think the plot we first shared is not great for illustrating it though, so apologies for that. The below plot should make our problem clearer (this is a different signal point/workspace to the above, and uses calc_kwargs = {'base_dir': 'clipped_normal'} ):

The CLs values, both expected and observed, jump up/down as we scan mu, rather than falling smoothly. Do you have any ideas what might be causing this?

[kratsg]

It's not clear why this behavior is happening in particular. Without access to the workspace, it's hard to really comment more.