docs: documentation of covariance and correlated fits extended.

This commit is contained in:
Fabian Joswig 2022-03-05 08:13:24 +00:00
parent a7ff26ed9c
commit 6bd3868179
3 changed files with 22 additions and 14 deletions

View file

@ -353,6 +353,11 @@ def func(a, x):
return a[0] * x1 ** 2 + a[1] * x2
```
`pyerrors` also supports correlated fits which can be triggered via the parameter `correlated_fit=True`.
Details about how the required covariance matrix is estimated can be found in `pyerrors.obs.covariance`.
Direct visualizations of the performed fits can be triggered via `resplot=True` or `qqplot=True`. For all available options see `pyerrors.fits.least_squares`.
## Total least squares fits
`pyerrors` can also fit data with errors on both the dependent and independent variables using the total least squares method also referred to orthogonal distance regression as implemented in [scipy](https://docs.scipy.org/doc/scipy/reference/odr.html), see `pyerrors.fits.least_squares`. The syntax is identical to the standard least squares case, the only diffrence being that `x` also has to be a `list` or `numpy.array` of `Obs`.

View file

@ -86,7 +86,7 @@ def least_squares(x, y, func, priors=None, silent=False, **kwargs):
```
It is important that all numpy functions refer to autograd.numpy, otherwise the differentiation
will not work
will not work.
priors : list, optional
priors has to be a list with an entry for every parameter in the fit. The entries can either be
Obs (e.g. results from a previous fit) or strings containing a value and an error formatted like
@ -95,24 +95,25 @@ def least_squares(x, y, func, priors=None, silent=False, **kwargs):
If true all output to the console is omitted (default False).
initial_guess : list
can provide an initial guess for the input parameters. Relevant for
non-linear fits with many parameters.
non-linear fits with many parameters.
method : str, optional
can be used to choose an alternative method for the minimization of chisquare.
The possible methods are the ones which can be used for scipy.optimize.minimize and
migrad of iminuit. If no method is specified, Levenberg-Marquard is used.
Reliable alternatives are migrad, Powell and Nelder-Mead.
resplot : bool
If true, a plot which displays fit, data and residuals is generated (default False).
qqplot : bool
If true, a quantile-quantile plot of the fit result is generated (default False).
expected_chisquare : bool
If true prints the expected chisquare which is
corrected by effects caused by correlated input data.
This can take a while as the full correlation matrix
has to be calculated (default False).
correlated_fit : bool
If true, use the full correlation matrix in the definition of the chisquare
(only works for prior==None and when no method is given, at the moment).
If True, use the full inverse covariance matrix in the definition of the chisquare cost function.
For details about how the covariance matrix is estimated see `pyerrors.obs.covariance`.
In practice the correlation matrix is Cholesky decomposed and inverted (instead of the covariance matrix).
This procedure should be numerically more stable as the correlation matrix is typically better conditioned (Jacobi preconditioning).
At the moment this option only works for `prior==None` and when no `method` is given.
expected_chisquare : bool
If True estimates the expected chisquare which is
corrected by effects caused by correlated input data (default False).
resplot : bool
If True, a plot which displays fit, data and residuals is generated (default False).
qqplot : bool
If True, a quantile-quantile plot of the fit result is generated (default False).
'''
if priors is not None:
return _prior_fit(x, y, func, priors, silent=silent, **kwargs)

View file

@ -1348,7 +1348,9 @@ def covariance(obs, visualize=False, correlation=False, **kwargs):
Notes
-----
The covariance is estimated by calculating the correlation matrix assuming no autocorrelation and then rescaling the correlation matrix by the full errors including the previous gamma method estimate for the autocorrelation of the observables. For observables defined on a single ensemble this is equivalent to assuming that the integrated autocorrelation time of an off-diagonal element is equal to the geometric mean of the integrated autocorrelation times of the corresponding diagonal elements.
The covariance is estimated by calculating the correlation matrix assuming no autocorrelation and then rescaling the correlation matrix by the full errors including the previous gamma method estimate for the autocorrelation of the observables. The covariance at windowsize 0 is guaranteed to be positive semi-definite
$$v_i\Gamma_{ij}(0)v_j=\frac{1}{N}\sum_{s=1}^N\sum_{i,j}v_i\delta_i^s\delta_j^s v_j=\frac{1}{N}\sum_{s=1}^N\sum_{i}|v_i\delta_i^s|^2\geq 0\,,$$ for every $v_i\in\mathbb{R}^N$, while such an identity does not hold for larger windows/lags.
For observables defined on a single ensemble our approximation is equivalent to assuming that the integrated autocorrelation time of an off-diagonal element is equal to the geometric mean of the integrated autocorrelation times of the corresponding diagonal elements.
$$\tau_{\mathrm{int}, ij}=\sqrt{\tau_{\mathrm{int}, i}\times \tau_{\mathrm{int}, j}}$$
This construction ensures that the estimated covariance matrix is positive semi-definite (up to numerical rounding errors).
'''