diff --git a/pyerrors/fits.py b/pyerrors/fits.py index 6518b224..b58f62b6 100644 --- a/pyerrors/fits.py +++ b/pyerrors/fits.py @@ -365,7 +365,7 @@ def least_squares(x, y, func, priors=None, silent=False, **kwargs): raise Exception('The minimization procedure did not converge.') output.chisquare = chisquare - output.dof = x_all.shape[-1] - n_parms + len(loc_priors) + output.dof = y_all.shape[-1] - n_parms + len(loc_priors) output.p_value = 1 - scipy.stats.chi2.cdf(output.chisquare, output.dof) if output.dof > 0: output.chisquare_by_dof = output.chisquare / output.dof @@ -393,7 +393,7 @@ def least_squares(x, y, func, priors=None, silent=False, **kwargs): hat_vector = prepare_hat_matrix() A = W @ hat_vector P_phi = A @ np.linalg.pinv(A.T @ A) @ A.T - expected_chisquare = np.trace((np.identity(x_all.shape[-1]) - P_phi) @ W @ cov @ W) + expected_chisquare = np.trace((np.identity(y_all.shape[-1]) - P_phi) @ W @ cov @ W) output.chisquare_by_expected_chisquare = output.chisquare / expected_chisquare if not silent: print('chisquare/expected_chisquare:', output.chisquare_by_expected_chisquare) diff --git a/tests/fits_test.py b/tests/fits_test.py index dbc227dc..3ed8f5fb 100644 --- a/tests/fits_test.py +++ b/tests/fits_test.py @@ -1142,6 +1142,53 @@ def test_fit_dof(): assert cd[0] != cd[0] # Check for nan assert np.all(np.array(cd[1:]) > 0) + N = 5 + + def fitf(a, x): + return a[0] + 0 * x + + def fitf_multi(a, x): + return a[0] + 0 * x[0] + 0*x[1] + + for priors in [None, [pe.cov_Obs(3, 1, 'p')]]: + if priors is None: + lp = 0 + else: + lp = len(priors) + x = [1. for i in range(N)] + y = [pe.cov_Obs(i, .1, '%d' % (i)) for i in range(N)] + [o.gm() for o in y] + res = pe.fits.least_squares(x, y, fitf, expected_chisquare=True, priors=priors) + assert(res.dof == N - 1 + lp) + if priors is None: + assert(np.isclose(res.chisquare_by_expected_chisquare, res.chisquare_by_dof)) + + kl = ['a', 'b'] + x = {k: [1. for i in range(N)] for k in kl} + y = {k: [pe.cov_Obs(i, .1, '%d%s' % (i, k)) for i in range(N)] for k in kl} + [[o.gm() for o in y[k]] for k in y] + res = pe.fits.least_squares(x, y, {k: fitf for k in kl}, expected_chisquare=True, priors=priors) + assert(res.dof == 2 * N - 1 + lp) + if priors is None: + assert(np.isclose(res.chisquare_by_expected_chisquare, res.chisquare_by_dof)) + + x = np.array([[1., 2.] for i in range(N)]).T + y = [pe.cov_Obs(i, .1, '%d' % (i)) for i in range(N)] + [o.gm() for o in y] + res = pe.fits.least_squares(x, y, fitf_multi, expected_chisquare=True, priors=priors) + assert(res.dof == N - 1 + lp) + if priors is None: + assert(np.isclose(res.chisquare_by_expected_chisquare, res.chisquare_by_dof)) + + x = {k: np.array([[1., 2.] for i in range(N)]).T for k in kl} + y = {k: [pe.cov_Obs(i, .1, '%d%s' % (i, k)) for i in range(N)] for k in kl} + [[o.gm() for o in y[k]] for k in y] + res = pe.fits.least_squares(x, y, {k: fitf_multi for k in kl}, expected_chisquare=True, priors=priors) + + assert(res.dof == 2 * N - 1 + lp) + if priors is None: + assert(np.isclose(res.chisquare_by_expected_chisquare, res.chisquare_by_dof)) + def test_combined_fit_constant_shape(): N1 = 16