lattice/pyerrors
1
0
Fork 0
mirror of https://github.com/fjosw/pyerrors.git synced 2025-03-15 06:40:24 +01:00
Code Issues Projects Releases Packages Wiki Activity

Merge pull request #85 from fjosw/fix/pinv_for_chi_exp

Browse source 
Replaced inverse by the pseudo inverse in the calculation of chi_exp
This commit is contained in:
Fabian Joswig 2022-03-07 11:59:31 +00:00 • committed by GitHub
commit e5d7271a2b
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
2 changed files with 19 additions and 2 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

4
pyerrors/fits.py
View file

@ -252,7 +252,7 @@ def total_least_squares(x, y, func, silent=False, **kwargs):
new_jac = np.concatenate((fused_row1, fused_row2), axis=1)
A = W @ new_jac
P_phi = A @ np.linalg.inv(A.T @ A) @ A.T
P_phi = A @ np.linalg.pinv(A.T @ A) @ A.T
expected_chisquare = np.trace((np.identity(P_phi.shape[0]) - P_phi) @ W @ cov @ W)
if expected_chisquare <= 0.0:
warnings.warn("Negative expected_chisquare.", RuntimeWarning)
@ -531,7 +531,7 @@ def _standard_fit(x, y, func, silent=False, **kwargs):
W = np.diag(1 / np.asarray(dy_f))
cov = covariance(y)
A = W @ jacobian(func)(fit_result.x, x)
P_phi = A @ np.linalg.inv(A.T @ A) @ A.T
P_phi = A @ np.linalg.pinv(A.T @ A) @ A.T
expected_chisquare = np.trace((np.identity(x.shape[-1]) - P_phi) @ W @ cov @ W)
output.chisquare_by_expected_chisquare = chisquare / expected_chisquare
if not silent:

17
tests/fits_test.py
View file

@ -319,6 +319,23 @@ def test_prior_fit():
fitp = pe.fits.least_squares([0, 1], y, f, priors=y, resplot=True, qqplot=True)
def test_correlated_fit_covobs():
x1 = pe.cov_Obs(1.01, 0.01 ** 2, 'test1')
x2 = pe.cov_Obs(2.01, 0.01 ** 2, 'test2')
x3 = pe.cov_Obs(2.99, 0.01 ** 2, 'test3')
[o.gamma_method() for o in [x1, x2, x3]]
def func(a, x):
return a[0] * x + a[1]
fit_res = pe.fits.least_squares(np.arange(1, 4), [x1, x2, x3], func, expected_chisquare=True)
assert np.isclose(fit_res.chisquare_by_dof, fit_res.chisquare_by_expected_chisquare)
fit_res_corr = pe.fits.least_squares(np.arange(1, 4), [x1, x2, x3], func, correlated_fit=True)
assert np.isclose(fit_res.chisquare_by_dof, fit_res_corr.chisquare_by_dof)
def test_error_band():
def f(a, x):
return a[0] + a[1] * x