From d5d42f184194052444c160b11609c266b9c72745 Mon Sep 17 00:00:00 2001 From: William Setterberg Date: Tue, 26 Nov 2024 13:20:46 -0600 Subject: [PATCH] fit range on levenberg-marquadt --- yaff/fitting.py | 17 +++++++++++++++-- 1 file changed, 15 insertions(+), 2 deletions(-) diff --git a/yaff/fitting.py b/yaff/fitting.py index 451d632..2723f32 100644 --- a/yaff/fitting.py +++ b/yaff/fitting.py @@ -538,11 +538,20 @@ def free_param_units(self) -> list[u.Unit]: return list(p.unit for p in self.free_parameters) -def levenberg_minimize(fitter: BayesFitter, **scipy_kwargs) -> BayesFitter: +def levenberg_minimize( + fitter: BayesFitter, + restriction: np.ndarray[bool] = None, + **scipy_kwargs +) -> BayesFitter: """Given a Bayes fitter, minimize its parameters using the Levenberg-Marquadt (weighted) least squares minimization like XSPEC does. - This minimization technique operates on **all** of the model and data + The `restriction` parameter is used to optionally exclude count bins + which shouldn't be considered when fitting. + + --- + + Levenberg-Marquadt operates on **all** of the model and data count bins. So, it tends to be more robust (and converge faster) than algorithms which are based on a single summary "fit statistic." @@ -550,6 +559,9 @@ def levenberg_minimize(fitter: BayesFitter, **scipy_kwargs) -> BayesFitter: as it will (robustly) put the parameter vector near its global minimum. """ + if restriction is None: + restriction = np.ones_like(fitter.data.counts, dtype=bool) + def residual_function(vector: ArrayLike): fitter.emplace_free_parameters(vector) if np.isneginf(fitter.eval_priors()): @@ -566,6 +578,7 @@ def residual_function(vector: ArrayLike): ) ret = (mod - compare) / total_error + ret[~restriction] = 0 # Any "zero-error" bins need to get deleted return np.nan_to_num(ret, copy=False, nan=0, posinf=0, neginf=0)