ENH: add optional argument for scaled covariance matrix of ODR

[AtsushiSakai]

Reference issue

fix #6842

What does this implement/fix?

fix covariance matrix of ODR to be scaled by the residual variance and added its tests.
I added optional argument for scaled covariance matrix of ODR.
I updated sd_beta docstring to clarify the value is scaled.

[pv]

See discussion in https://mail.python.org/pipermail/scipy-user/2013-February/034196.html and #448 on probably similar resolved issue in curve_fit --- there it appears people from different fields have conflicting expectations on what's the right thing to do with the covariance matrix. If this issue is the same as there, the solution probably should be the same. So what should be done first is to check whether this indeed is the same issue.

[AtsushiSakai]

Thank you so much for your comment. I didn’t realize that there are some different expectations for covariance matrix. I agree we should not break backward compatibility, there should be an optional argument to scale it or not. I will try to fix it.

[kurisutsukato]

I made a comment related to this issue here:
#6842 (comment)

Basically there seems to be an inconsitency between curve_fit and odr in mode 2 (=least squares) in terms of the calculation of the error of the parameters. The unweighted covariance matrix returned by odr is almost equal to that returned by curve_fit.

[AtsushiSakai]

@kurisutsukato Hi. Thank you for your comment. I commented in #6842. PTAL.

[juhuebner]

💡 @AtsushiSakai
I suggest to make the additional options and the related output comparable, or better: identical to scipy.optimize.curve_fit() . I followed the discussion here #6842 as well. There is, however, a possible issue with scipy.optimize.curve_fit() as well, since omitting the error data set sigma and putting absolute_sigma=True still gives an output, but does not make sense since because in curve_fit(func, x_data, y_data, absolute_sigma=True) the weighting, i.e., 1/σ is defaulted to unity and the option absolute_sigma=True makes scipy.optimize.curve_fit() think the real errors are all ones. The correct output would be only given by curve_fit(func, x_data, y_data) - where absolute_sigma is defaulted to False, i.e. the scaling of the covariance matrix matrix in this case, as it discussed here makes sense and one has to be careful what is scaled in which scenario.

see issue here: #12584