Sigma keyword for hesse

[HDembinski]

This addresses issue #173, by making hesse(...) and minos(...) accept a sigma keyword. The effect could be obtained in other ways, but most people do not understand how set_errordef(...) works, and I don't blame them, it is not obvious. I like the sigma keyword in minos(...), because it is pretty clear.

People expect similar things to behave similarly, and hesse and minos are pretty similar. It therefore makes sense to have a similar interface.

[cdeil]

@HDembinski - What effect does this have? I.e. which properties / results does this change influence?
Presumably m.errors? And will m.covariance change or not? Others?

I think the question here is whether this addition makes the overall situation more or less confusing. Probably an argument could be made either way. On the one hand, it's convenient to set sigma when calling hesse, especially given that minos already has a sigma. On the other hand, this introduces a second (or third?) way to do the same thing, i.e. could lead do more different usage and possibly more confusion among users.

I think whether we merge this addition or not, this is one of the points people (including me) struggle with most, and IMO the most important improvement would be docs additions, e.g. a page or section that explains this (or link to the one in Fred James writeup if available), and then to link to that explanation from all methods / properties that depend on this choice or "up" or "errordef" or "sigma".

@HDembinski, @piti118, or @jetteodder from #173 or anyone - Do you think we should add this? Any other thoughts on how to make this point easier / clearer for users?

[HDembinski]

@cdeil I value consistency very much. Yes, hesse(sigma=2) gives you consistently 2 sigma errors when you call errors() and covariance(). I tested that on my computer, but you are right I need to add a unit test for that. Before I polish this PR, I would like to get some feedback though, because of the design points you raised.

I think this interface change is a good idea. I studied the Zen of Python carefully. I think it is great design advice and one of the reasons why I love Python.
https://www.python.org/dev/peps/pep-0020/

Let's go through it and see what speaks for adding the sigma keyword.
"Explicit is better than implicit."
"Simple is better than complex."
"Readability counts."
Firstly, hesse(sigma=1) is more explicit (even if it just says so in the function signature) than hesse(), because the keyword in the signature says what kind of errors you get (1 sigma errors, not 2, 3 ...). This is self-documenting.
Secondly, hesse(sigma=2) is way simpler than

some_sigma = 2
old_up = m.get_errordef()
m.set_errordef(old_up * some_sigma**2)
m.hesse()
m.set_errordef(old_up)

This is how it is internally implemented for hesse and minos. For hesse but not for minos, we could even further optimize this internally by skipping the explicit re-calculation of hesse and just scale errors and covariance matrix.

"There should be one-- and preferably only one --obvious way to do it."
Someone already decided in the past that a sigma keyword is nice to have for minos(). I think the decision to add it to minos() was a good one. It now follows out of consistency that hesse() should have it as well. Users expect consistency, and issue #173 proves this.

Consistency is so important in design, because it allows you to use the software without reading the documentation. Nobody likes to read the documentation. If an interface is consistent, you rarely have to go back to the documentation to figure out what is going on. It is easier to remember how to do things, because patterns emerge that you learn intuitively. Again, consistency is really really important. It is one of the key features of Python that people fall in love with. We want them to fall in love with iminuit, too, don't we? Btw, it is also one of the reasons why I hate using ROOT, because it is not consistent at all.

"Now is better than never."
"Although never is often better than right now."
We need to discuss interface changes such as this one, I absolutely agree. But eventually we should make a decision and move on. The two proposed solutions to #173 are

1. make the interface consistent
2. keep the interface unchanged (and inconsistent) and document the issue

[cdeil]

@HDembinski - I don't have an opinion on this change yet, I'm about 50:50 on the proposal. I just want to be careful with changes in iminuit, and if we make them, try to get them as "right" as possible, i.e. do as little changes over the years as possible (given that Minuit is 40? years old and iminuit hasn't been changed in a few years).

---

One concern with adding sigma to hesse is that people will call it multiple times to get 1, 2, 3 sigma errors. That is very computationally inefficient, no? Wouldn't it be better to teach via the docs how to get 1, 2, 3 sigma errors, just by multiplying m.errors times sigma? In that sense, hesse is not equivalent to minos, where you have to re-run to get a different sigma.

Another concern is that many people just call migrad and after that access errors, even if it's not recommended. So if consistency is an argument, you'd have to add sigma=1 to migrad as well.

My strongest concern here is about m.covariance. Is it affected by the sigma you pass to hesse? It shouldn't be, right? It should always be the same, and from the one covariance matrix one can compute errors at different sigma levels, no?

Just for reference, this is what we have for errordef docs currently:
http://iminuit.readthedocs.io/en/latest/api.html?highlight=errordef#iminuit.Minuit
http://iminuit.readthedocs.io/en/latest/api.html?highlight=errordef#iminuit.Minuit.set_up

Given that already a long time ago this point was mentioned as the number 1 error people make when fitting with Minuit (or in general),
http://lmu.web.psi.ch/docu/manuals/software_manuals/minuit2/mnerror.pdf
our docs in iminuit are falling short.

---

@HDembinski - Please comment on the concerns I mentioned, especially the one about whether the sigma to HESSE affects covariance or not, and if you think it should or not. (also, the question if this should then also be added to migrad or not given that it estimates and returns an error is important to have a complete proposal on the table)

@piti118 and all - please give your opinion here.

[HDembinski]

@cdeil Please find below answers to all the points you raised.

I just want to be careful with changes in iminuit, and if we make them, try to get them as "right" as possible, i.e. do as little changes over the years as possible (given that Minuit is 40? years old and iminuit hasn't been changed in a few years).

This is not a change in Minuit itself. This is a useful convenience feature of the iminuit Wrapper. Neither the original minos(...) nor the original hesse(...) functions in C++ have a sigma parameter. This is something that was added in the iminuit wrapper.

One concern with adding sigma to hesse is that people will call it multiple times to get 1, 2, 3 sigma errors. That is very computationally inefficient, no? Wouldn't it be better to teach via the docs how to get 1, 2, 3 sigma errors, just by multiplying m.errors times sigma? In that sense, hesse is not equivalent to minos, where you have to re-run to get a different sigma.

No, it can be made as efficient as doing it by hand. Let me quote myself:

For hesse but not for minos, we could even further optimize this internally by skipping the explicit re-calculation of hesse and just scale errors and covariance matrix.

Also, performance is secondary to usability in Python. If you want maximum performance, you wouldn't use Python in the first place.

My strongest concern here is about m.covariance. Is it affected by the sigma you pass to hesse?

Yes, it is, I already said so above.

It shouldn't be, right? It should always be the same, and from the one covariance matrix one can compute errors at different sigma levels, no?

It should be affected. The normal covariance matrix is what you use to error propagate 1 sigma errors. The other covariance matrices can be used to propagate 2 or 3 or N sigma errors. People who want to use this feature will know what to do with it. Everyone else will just use the default, which is sigma=1 and never bother.

Just for reference, this is what we have for errordef docs currently:
http://iminuit.readthedocs.io/en/latest/api.html?highlight=errordef#iminuit.Minuit

The documentation is bad. The true purpose of errordef in __init__ is to choose between the values 1 for least-squares functions and 0.5 for likelihood functions. I will address this problem in another PR.

http://iminuit.readthedocs.io/en/latest/api.html?highlight=errordef#iminuit.Minuit.set_up

This is just a forward, which I do not like at all. set_up is even more cryptic than set_errordef, and it is just a forward to set_errordef. The documentation in set_errordef correctly refers to the purpose that errordef has in __init__, the explanation needs to be copied there. That would be a separate PR.

Given that already a long time ago this point was mentioned as the number 1 error people make when fitting with Minuit (or in general),
http://lmu.web.psi.ch/docu/manuals/software_manuals/minuit2/mnerror.pdf
our docs in iminuit are falling short.

This I can easily believe without even knowing any statistics about how often users make certain errors.

[HDembinski]

As a side remark, I know that trust in someone requires some experience and time, but to raise my credibility somewhat, I would like to say that I spend a good amount of my free time thinking about interface design, and I am very aware of the long-term problems caused by bad design decisions.

I am the author of the histogram library, which will be reviewed by the Boost community in a few weeks to decide about including it in Boost. Once the library is in Boost, the interface is fixed and cannot be changed, so I am thinking constantly about how to not make decisions which are limiting us later on.

I read several books about design patterns in C++ from Stroustrup, Sutter, Alexandrescu, and Meyers, and Gamma. C++ is a different language, but conceptually the challenges of good design are the same in both languages, Python just uses different idioms.

[HDembinski]

@piti118 What do you think about this change? @cdeil and I cannot agree on our own whether this is a good idea or not, we would be happy to get your opinion. :)

[cdeil]

@HDembinski and I discussed this at lunch yesterday, and now that I understand the implications better, here's what I think:

I'm -1 on adding sigma to hesse.

I admit that it can be a nice convenience if someone wants "2 sigma errors" in their script or interactive analysis session. But I think that when collaborating on scripts with others or re-using scripts from others, this addition will cause problems. Currently Minuit.errors and Minuit.covariance are always "1 sigma" errors. With this "sigma" option in "hesse", "errors" and "covariance" will be "1 sigma" errors in many cases (currently all), but "n sigma" errors sometimes if users start to pass "sigma=n" to "hesse". I find this confusing, especially for the covariance matrix. To also quote from the zen of Python: "There should be one-- and preferably only one --obvious way to do it."
For iminuit IMO that's what we currently have: the user runs "hesse" and accesses "errors" and "covariance", and if they want 2, 3, 4 sigma errors, they multiply the errors with 2, 3, 4 after and the covariance with 4, 9, 16. That's also more efficient than multiple calls into hesse, and more consistent with the fact that "migrad" also estimates the "errors" and "covariance" properties and sets them to 1-sigma always. You're proposing to add a second equivalent way to do the same thing. IMO the addition of sigma in minos was justified, as there to get 1 and 2 sigma errors, the minos algorithm needs to run twice. That's not the case with the hesse algorithm.

@HDembinski - I think the drawbacks of this addition would become more obvious if the PR was complete in the sense that you'd try to accurately describe what is returned in the docstrings of Minuit.errors and Minuit.covariance. It would have to be changed to something like "usually these are "1 sigma errors" or "the normal covariance", but if you called "hesse(sigma=n)", then it's scaled by "n" (for error) or "n^2" (for covariance) already.", no?

Just to be clear: I think this proposal and discussion was very useful. I'm just -1 on improving iminuit as proposed here, and the counter-proposal is to improve the docs how to scale errors and covariance instead (and to put in the numpy array versions of those, to make it even easier, as discussed in another issue).

[HDembinski]

@cdeil I can counter all your arguments.

1. You are mixing two questions, one is about the interface and one is about performance. Multiple calls to hesse with sigma=N can be made as performant like a manual handling, if we internally track these calls. This is easy to do. Therefore, this is not a show-stopper. Let's not confuse questions about the interface with the implementation.

2. You worry about inconsistency of the meaning of errors and covariance when hesse is called with keyword sigma=2. We already have this situation in merrors. This attribute has a different meaning already, depending on whether you call minos with sigma=1 or sigma=2. Whoever is actually calling hesse or minos with sigma=2 knows what they are doing. For anyone else the attributes have the expected meaning. There is no problem here, and this change is 100 % backward compatible.

3. You have attacked none of my arguments concerning interface consistency and the new self-documenting property of the sigma=1 keyword in the hesse signature. Therefore, these arguments still stand.

When you form an opinion, please fairly consider all the arguments, all pros and cons. That's what I do. I have countered your cons, and there are still pros left. So the rational choice is to accept this change.

[cdeil]

@HDembinski - if you really still think this is a good proposal and want it in, please make the pull request complete so that the implications become more clear. At the very least you need to amend the docstring of the errors and covariance properties to explain what they will now contain.

Multiple calls to hesse with sigma=N can be made as performant like a manual handling, if we internally track these calls. This is easy to do.

I don't think we should add such complexity (more ways to do the same thing and state tracking than we already have, i.e. triggering hesse or migrad from values or errors property access). If this is your vision for what to do as long-term solution, I would suggest you do that change here, to show what it looks like. I think it's very relevant to this proposal / decision.

As it stands, I think this PR and advertising to users to do multiple calls to hesse to get 1, 2, 3 sigma errors is terrible advice, if they can just do this:

errors_1_sigma = 1 * m.np_errors()
errors_2_sigma = 2 * m.np_errors()
errors_3_sigma = 3 * m.np_errors()

@piti118 - Please comment.

[HDembinski]

Closing this request, because @cdeil and I could not agree whether this is a useful feature or not. We will add documentation instead to show how to handle this.