diff --git a/src/easyscience/fitting/engine_base.py b/src/easyscience/fitting/engine_base.py new file mode 100644 index 00000000..08190ea8 --- /dev/null +++ b/src/easyscience/fitting/engine_base.py @@ -0,0 +1,200 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause + +from abc import ABCMeta +from inspect import Parameter as InspectParameter +from inspect import Signature +from inspect import _empty +from typing import Callable +from typing import Dict +from typing import Tuple + +import numpy as np + +# causes circular import when Parameter is imported +# from easyscience.base_classes import ObjBase +from easyscience.variable import Parameter + +PARAMETER_PREFIX = 'p' + + +class EngineBase(metaclass=ABCMeta): + """ + Base for all evaluation engines — minimizers and samplers alike. + + An engine binds an EasyScience object and a fit function, and + repeatedly evaluates the function while writing values back into the + object's ``Parameter`` instances. ``EngineBase`` owns that shared + machinery: the parameter cache, the ``Parameter``-writing wrapped + fit function, and value restore on failure. It deliberately declares + no abstract methods: the run interfaces live on its subclasses + (``MinimizerBase.fit``, ``DreamSampler.run``). + """ + + package: str = None + + def __init__( + self, + obj, #: ObjBase, + fit_function: Callable, + ): # todo after constraint changes, add type hint: obj: ObjBase # noqa: E501 + self._object = obj + self._original_fit_function = fit_function + self._cached_pars: Dict[str, Parameter] = {} + self._cached_pars_vals: Dict[str, Tuple[float, float]] = {} + self._fit_function = None + + def _restore_parameter_values(self) -> None: + for key in self._cached_pars.keys(): + self._cached_pars[key].value = self._cached_pars_vals[key][0] + self._cached_pars[key].error = self._cached_pars_vals[key][1] + + def evaluate( + self, x: np.ndarray, minimizer_parameters: dict[str, float] | None = None, **kwargs + ) -> np.ndarray: + """ + Evaluate the fit function for values of x. + + Parameters used are either the latest or user supplied. If the + parameters are user supplied, it must be in a dictionary of + {'parameter_name': parameter_value,...}. + + Parameters + ---------- + x : np.ndarray + X values for which the fit function will be evaluated. + minimizer_parameters : dict[str, float] | None, default=None + Dictionary of parameters which will be used in the fit + function. They must be in a dictionary of {'parameter_name': + parameter_value,...}. By default, None. + **kwargs : + Additional arguments. + + Returns + ------- + np.ndarray + Y values calculated at points x for a set of parameters. + + Raises + ------ + TypeError + If ``minimizer_parameters`` is not a dictionary. + """ + if minimizer_parameters is None: + minimizer_parameters = {} + if not isinstance(minimizer_parameters, dict): + raise TypeError('minimizer_parameters must be a dictionary') + + if self._fit_function is None: + # This will also generate self._cached_pars + self._fit_function = self._generate_fit_function() + + minimizer_parameters = self._prepare_parameters(minimizer_parameters) + + return self._fit_function(x, **minimizer_parameters, **kwargs) + + def _prepare_parameters(self, parameters: dict[str, float]) -> dict[str, float]: + """ + Prepare the parameters for the engine. + + Parameters + ---------- + parameters : dict[str, float] + Dict of parameters for the engine with names as keys. + + Returns + ------- + dict[str, float] + Completed parameter dictionary for the engine. + """ + pars = self._cached_pars + + for name, item in pars.items(): + parameter_name = PARAMETER_PREFIX + str(name) + if parameter_name not in parameters.keys(): + parameters[parameter_name] = item.value + return parameters + + def _generate_fit_function(self) -> Callable: + """ + Using the user supplied ``fit_function``, wrap it in such a way + we can update ``Parameter`` on iterations. + + Returns + ------- + Callable + A fit function which is compatible with bumps models. + """ + # Original fit function + func = self._original_fit_function + # Get a list of `Parameters` + self._cached_pars = {} + self._cached_pars_vals = {} + for parameter in self._object.get_fit_parameters(): + key = parameter.unique_name + self._cached_pars[key] = parameter + self._cached_pars_vals[key] = (parameter.value, parameter.error) + + # Make a new fit function + def _fit_function(x: np.ndarray, **kwargs) -> np.ndarray: + """ + Wrapped fit function which now has an EasyScience compatible + form. + + Parameters + ---------- + x : np.ndarray + Array of data points to be calculated. + **kwargs : + Key word arguments. + + Returns + ------- + np.ndarray + Points calculated at ``x``. + """ + # Update the `Parameter` values and the callback if needed + # TODO THIS IS NOT THREAD SAFE :-( + + for name, value in kwargs.items(): + par_name = name[1:] + if par_name in self._cached_pars.keys(): + # This will take into account constraints + if self._cached_pars[par_name].value != value: + self._cached_pars[par_name].value = value + + # Since we are calling the parameter fset will be called. + # TODO Pre processing here + return_data = func(x) + # TODO Loading or manipulating data here + return return_data + + _fit_function.__signature__ = self._create_signature(self._cached_pars) + return _fit_function + + @staticmethod + def _create_signature(parameters: Dict[int, Parameter]) -> Signature: + """ + Wrap the function signature. + + This is done as lmfit wants the function to be in the form: f = + (x, a=1, b=2)... Where we need to be generic. Note that this + won't hold for much outside of this scope. + """ + wrapped_parameters = [] + wrapped_parameters.append( + InspectParameter('x', InspectParameter.POSITIONAL_OR_KEYWORD, annotation=_empty) + ) + + for name, parameter in parameters.items(): + default_value = parameter.value + + wrapped_parameters.append( + InspectParameter( + PARAMETER_PREFIX + str(name), + InspectParameter.POSITIONAL_OR_KEYWORD, + annotation=_empty, + default=default_value, + ) + ) + return Signature(wrapped_parameters) diff --git a/src/easyscience/fitting/fitter.py b/src/easyscience/fitting/fitter.py index 7cb948e0..e8e0d4d8 100644 --- a/src/easyscience/fitting/fitter.py +++ b/src/easyscience/fitting/fitter.py @@ -283,7 +283,7 @@ def inner_fit_callable( y: np.ndarray, weights: Optional[np.ndarray] = None, vectorized: bool = False, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, **kwargs, ) -> FitResults: """ @@ -418,112 +418,3 @@ def _post_compute_reshaping( fit_result.y_calc = np.reshape(fit_result.y_calc, y.shape) fit_result.y_err = np.reshape(fit_result.y_err, y.shape) return fit_result - - def mcmc_sample( - self, - x: np.ndarray, - y: np.ndarray, - weights: np.ndarray, - samples: int = 10000, - burn: int = 2000, - thin: int = 10, - population: Optional[int] = None, - vectorized: bool = False, - sampler_kwargs: Optional[dict] = None, - progress_callback: Optional[Callable[[dict], Optional[bool]]] = None, - abort_test: Optional[Callable[[], bool]] = None, - ) -> dict: - """ - Run Bayesian MCMC sampling using the BUMPS DREAM sampler. - - Works with both a plain ``Fitter`` (single dataset) and a - ``MultiFitter`` (multiple datasets) via polymorphic dispatch: - ``_precompute_reshaping`` and ``_fit_function_wrapper`` are - resolved on the concrete subclass at call time, so multi-dataset - flattening is handled automatically when called on a - ``MultiFitter`` instance. - - Parameters - ---------- - x : np.ndarray - Independent variable array (or list of arrays for - ``MultiFitter``). - y : np.ndarray - Dependent variable array (or list of arrays for - ``MultiFitter``). - weights : np.ndarray - Weight array (or list of arrays for ``MultiFitter``). - samples : int, default=10000 - Number of retained DREAM samples requested from BUMPS. - burn : int, default=2000 - Burn-in steps to discard before collecting samples. - thin : int, default=10 - Thinning interval — only every ``thin``-th sample is kept, - which reduces autocorrelation between consecutive draws. - population : Optional[int], default=None - BUMPS DREAM population count (number of parallel chains). - vectorized : bool, default=False - When ``True``, each x array may be multi-dimensional (e.g. - an ``(N, M, 2)`` grid for a 2D model) and is left as-is. - When ``False`` (default), each x array is expected to be - 1-D. - sampler_kwargs : Optional[dict], default=None - Additional keyword arguments forwarded to the BUMPS DREAM - sampler. - progress_callback : Optional[Callable[[dict], Optional[bool]]], default=None - Optional callback invoked at each DREAM generation. The - payload dict includes ``iteration`` and ``sampling: True``. - abort_test : Optional[Callable[[], bool]], default=None - Optional callable that returns ``True`` to abort sampling - early. - - Returns - ------- - dict - Dictionary with keys ``'draws'``, ``'param_names'``, - ``'internal_bumps_object'``, and ``'logp'``. - - Raises - ------ - ValueError - If ``samples``, ``burn``, or ``thin`` are invalid. - RuntimeError - If the active minimizer is not a BUMPS instance. - """ - if not isinstance(samples, int) or samples <= 0: - raise ValueError('samples must be a positive integer.') - if not isinstance(burn, int) or burn < 0: - raise ValueError('burn must be a non-negative integer.') - if not isinstance(thin, int) or thin < 1: - raise ValueError('thin must be a positive integer.') - - x_fit, x_new, y_new, w_new, dims = self._precompute_reshaping(x, y, weights, vectorized) - self._dependent_dims = dims - - original_fit_func = self._fit_function - self.fit_function = self._fit_function_wrapper(x_new, flatten=True) - - try: - minimizer = self.minimizer - if not (hasattr(minimizer, 'package') and minimizer.package == 'bumps'): - raise RuntimeError( - 'Bayesian sampling requires a BUMPS minimizer. ' - 'Use ``fitter.switch_minimizer(AvailableMinimizers.Bumps)`` first.' - ) - - result = minimizer.mcmc_sample( - x=x_fit, - y=y_new, - weights=w_new, - samples=samples, - burn=burn, - thin=thin, - population=population, - sampler_kwargs=sampler_kwargs, - progress_callback=progress_callback, - abort_test=abort_test, - ) - finally: - self.fit_function = original_fit_func - - return result diff --git a/src/easyscience/fitting/minimizers/bumps_utils/__init__.py b/src/easyscience/fitting/minimizers/bumps_utils/__init__.py index 35255615..cdd3d6c6 100644 --- a/src/easyscience/fitting/minimizers/bumps_utils/__init__.py +++ b/src/easyscience/fitting/minimizers/bumps_utils/__init__.py @@ -2,6 +2,21 @@ # SPDX-License-Identifier: BSD-3-Clause from .eval_counter import EvalCounter +from .problem import build_curve_problem +from .problem import parameter_names +from .problem import parameter_snapshot +from .problem import to_bumps_parameter from .progress_monitor import BumpsProgressMonitor +from .validation import validate_arrays +from .validation import validate_run_settings -__all__ = ['BumpsProgressMonitor', 'EvalCounter'] +__all__ = [ + 'BumpsProgressMonitor', + 'EvalCounter', + 'build_curve_problem', + 'parameter_names', + 'parameter_snapshot', + 'to_bumps_parameter', + 'validate_arrays', + 'validate_run_settings', +] diff --git a/src/easyscience/fitting/minimizers/bumps_utils/problem.py b/src/easyscience/fitting/minimizers/bumps_utils/problem.py new file mode 100644 index 00000000..b5e25d42 --- /dev/null +++ b/src/easyscience/fitting/minimizers/bumps_utils/problem.py @@ -0,0 +1,142 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause +"""BUMPS problem construction shared by the ``Bumps`` minimizer and +``DreamSampler``. + +These are free functions rather than ``Bumps`` methods so that any +:class:`~easyscience.fitting.engine_base.EngineBase` — a minimizer or a +sampler — can build a BUMPS ``Curve``/``FitProblem`` without inheriting +from the minimizer. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING + +import numpy as np +from bumps.names import Curve +from bumps.names import FitProblem +from bumps.parameter import Parameter as BumpsParameter + +from easyscience.variable import Parameter + +from ...engine_base import PARAMETER_PREFIX +from .eval_counter import EvalCounter + +if TYPE_CHECKING: + from ...engine_base import EngineBase + + +def to_bumps_parameter(par: Parameter) -> BumpsParameter: + """Convert an EasyScience ``Parameter`` to a prefixed ``BumpsParameter``. + + Parameters + ---------- + par : Parameter + EasyScience parameter to convert. + + Returns + ------- + BumpsParameter + Bumps Parameter compatible object, named + ``PARAMETER_PREFIX + par.unique_name``. + """ + return BumpsParameter( + name=PARAMETER_PREFIX + par.unique_name, + value=par.value, + bounds=[par.min, par.max], + fixed=par.fixed, + ) + + +def build_curve_problem( + engine: 'EngineBase', + x: np.ndarray, + y: np.ndarray, + weights: np.ndarray, + parameters: list[Parameter] | None = None, +) -> tuple[FitProblem, EvalCounter, Curve]: + """Build a BUMPS ``FitProblem`` around an engine's wrapped fit function. + + Wraps ``engine._generate_fit_function()`` in an :class:`EvalCounter`, + converts the engine's cached parameters (or the explicitly supplied + ``parameters``) via :func:`to_bumps_parameter`, and assembles + ``Curve(fit_func, x, y, dy=1/weights, **bumps_pars)`` into a + ``FitProblem``. + + Parameters + ---------- + engine : EngineBase + The engine (minimizer or sampler) supplying the fit function and + parameter cache. + x : np.ndarray + Independent variable array. + y : np.ndarray + Dependent variable array. + weights : np.ndarray + Weight array; converted to ``dy = 1 / weights``. + parameters : list[Parameter] | None, default=None + Optional explicit EasyScience parameters to bind into the model + instead of the engine's cached parameters. + + Returns + ------- + tuple[FitProblem, EvalCounter, Curve] + The assembled problem, the evaluation counter wrapping the fit + function (exposes ``count`` for evaluation bookkeeping), and the + ``Curve`` model itself. The ``Curve`` is surfaced directly + because ``FitProblem.fitness`` is deprecated in BUMPS (>= 1.0.4 + it emits a ``UserWarning``) — callers must not go through it. + """ + fit_func = EvalCounter(engine._generate_fit_function()) + + bumps_pars = {} + if not parameters: + for name, par in engine._cached_pars.items(): + bumps_pars[PARAMETER_PREFIX + str(name)] = to_bumps_parameter(par) + else: + for par in parameters: + bumps_pars[PARAMETER_PREFIX + par.unique_name] = to_bumps_parameter(par) + + curve = Curve(fit_func, x, y, dy=1 / weights, **bumps_pars) + return FitProblem(curve), fit_func, curve + + +def parameter_names(problem: FitProblem) -> list[str]: + """Return the problem's parameter names with the prefix stripped. + + Parameters + ---------- + problem : FitProblem + A BUMPS problem built by :func:`build_curve_problem`. + + Returns + ------- + list[str] + Parameter names in problem order, without ``PARAMETER_PREFIX``. + """ + return [(p.name or '')[len(PARAMETER_PREFIX) :] for p in problem._parameters] + + +def parameter_snapshot(problem: FitProblem, point: np.ndarray | None) -> dict: + """Snapshot the problem's parameter values as ``{name: value}``. + + Parameters + ---------- + problem : FitProblem + A BUMPS problem built by :func:`build_curve_problem`. + point : np.ndarray | None + Parameter values to report; when ``None`` the problem's current + values (``problem.getp()``) are used. + + Returns + ------- + dict + Mapping of prefix-stripped parameter names to ``float`` values. + """ + labels = problem.labels() + values = problem.getp() if point is None else point + snapshot = {} + for label, value in zip(labels, values): + snapshot[label[len(PARAMETER_PREFIX) :]] = float(value) + return snapshot diff --git a/src/easyscience/fitting/minimizers/bumps_utils/validation.py b/src/easyscience/fitting/minimizers/bumps_utils/validation.py new file mode 100644 index 00000000..cb254a3b --- /dev/null +++ b/src/easyscience/fitting/minimizers/bumps_utils/validation.py @@ -0,0 +1,86 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause +"""Input validation shared by the BUMPS minimizer and the DREAM sampler.""" + +from __future__ import annotations + +import numpy as np + + +def validate_run_settings(samples: int, burn: int, thin: int) -> None: + """Validate the DREAM run settings. + + Parameters + ---------- + samples : int + Number of raw samples to draw; must be a positive integer. + burn : int + Burn-in generations to discard; must be a non-negative integer. + thin : int + Thinning interval; must be a positive integer. + + Raises + ------ + ValueError + If any value is out of range or not an integer. + """ + # bool is a subclass of int, so ``samples=True`` would otherwise pass as + # ``samples=1``; these checks are strict (``10.0`` is rejected), so + # booleans must be rejected too. + if not isinstance(samples, int) or isinstance(samples, bool) or samples <= 0: + raise ValueError('samples must be a positive integer.') + if not isinstance(burn, int) or isinstance(burn, bool) or burn < 0: + raise ValueError('burn must be a non-negative integer.') + if not isinstance(thin, int) or isinstance(thin, bool) or thin < 1: + raise ValueError('thin must be a positive integer.') + + +def validate_arrays( + x: np.ndarray, + y: np.ndarray, + weights: np.ndarray, + *, + check_finite_xy: bool = True, +) -> None: + """Validate the (x, y, weights) arrays for a BUMPS problem. + + Checks shape agreement between the three arrays, finiteness and + strict positivity of the weights, and — when ``check_finite_xy`` is + ``True`` — finiteness of x and y. Sampling passes ``True``; the + classical fit path passes ``False`` to keep its historically more + permissive behaviour. + + Parameters + ---------- + x : np.ndarray + Independent variable array. + y : np.ndarray + Dependent variable array. + weights : np.ndarray + Weight array (converted to ``dy = 1 / weights`` downstream). + check_finite_xy : bool, default=True + Also require x and y to be free of NaN/infinite values. + + Raises + ------ + ValueError + If the shapes disagree, the weights are non-finite or + non-positive, or (with ``check_finite_xy``) x/y are non-finite. + """ + if y.shape != x.shape: + raise ValueError('x and y must have the same shape.') + + if check_finite_xy: + if not np.isfinite(x).all(): + raise ValueError('x cannot contain NaN or infinite values.') + if not np.isfinite(y).all(): + raise ValueError('y cannot contain NaN or infinite values.') + + if weights.shape != x.shape: + raise ValueError('Weights must have the same shape as x and y.') + + if not np.isfinite(weights).all(): + raise ValueError('Weights cannot be NaN or infinite.') + + if (weights <= 0).any(): + raise ValueError('Weights must be strictly positive and non-zero.') diff --git a/src/easyscience/fitting/minimizers/minimizer_base.py b/src/easyscience/fitting/minimizers/minimizer_base.py index 48eef544..136e4396 100644 --- a/src/easyscience/fitting/minimizers/minimizer_base.py +++ b/src/easyscience/fitting/minimizers/minimizer_base.py @@ -1,17 +1,10 @@ # SPDX-FileCopyrightText: 2026 EasyScience contributors # SPDX-License-Identifier: BSD-3-Clause -from abc import ABCMeta from abc import abstractmethod -from inspect import Parameter as InspectParameter -from inspect import Signature -from inspect import _empty from typing import Any from typing import Callable -from typing import Dict from typing import List -from typing import Tuple -from typing import Union import numpy as np @@ -20,20 +13,23 @@ from easyscience.variable import Parameter from ..available_minimizers import AvailableMinimizers +from ..engine_base import PARAMETER_PREFIX +from ..engine_base import EngineBase from .utils import FitError from .utils import FitResults -MINIMIZER_PARAMETER_PREFIX = 'p' +# Back-compat alias: the canonical constant now lives in +# ``easyscience.fitting.engine_base`` and is shared by minimizers and +# samplers alike. +MINIMIZER_PARAMETER_PREFIX = PARAMETER_PREFIX -class MinimizerBase(metaclass=ABCMeta): +class MinimizerBase(EngineBase): """ This template class is the basis for all minimizer engines in ``EasyScience``. """ - package: str = None - def __init__( self, obj, #: ObjBase, @@ -42,14 +38,10 @@ def __init__( ): # todo after constraint changes, add type hint: obj: ObjBase # noqa: E501 if minimizer_enum.method not in self.supported_methods(): raise FitError(f'Method {minimizer_enum.method} not available in {self.__class__}') - self._object = obj - self._original_fit_function = fit_function + super().__init__(obj=obj, fit_function=fit_function) self._minimizer_enum = minimizer_enum self._method = minimizer_enum.method - self._cached_pars: Dict[str, Parameter] = {} - self._cached_pars_vals: Dict[str, Tuple[float]] = {} self._cached_model = None - self._fit_function = None @property def enum(self) -> AvailableMinimizers: @@ -59,11 +51,6 @@ def enum(self) -> AvailableMinimizers: def name(self) -> str: return self._minimizer_enum.name - def _restore_parameter_values(self) -> None: - for key in self._cached_pars.keys(): - self._cached_pars[key].value = self._cached_pars_vals[key][0] - self._cached_pars[key].error = self._cached_pars_vals[key][1] - @abstractmethod def fit( self, @@ -75,7 +62,7 @@ def fit( method: str | None = None, tolerance: float | None = None, max_evaluations: int | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, **kwargs, ) -> FitResults: """ @@ -99,8 +86,10 @@ def fit( Requested convergence tolerance. By default, None. max_evaluations : int | None, default=None Maximum number of objective evaluations. By default, None. - progress_callback : Callable[[dict], bool | None] | None, default=None - Optional progress callback. By default, None. + progress_callback : Callable[[dict], None] | None, default=None + Optional progress callback. Its return value is ignored by + every backend; use ``abort_test`` to stop a running fit. By + default, None. **kwargs : Additional arguments for the fitting function. @@ -110,50 +99,6 @@ def fit( Fit results. """ - def evaluate( - self, x: np.ndarray, minimizer_parameters: dict[str, float] | None = None, **kwargs - ) -> np.ndarray: - """ - Evaluate the fit function for values of x. - - Parameters used are either the latest or user supplied. If the - parameters are user supplied, it must be in a dictionary of - {'parameter_name': parameter_value,...}. - - Parameters - ---------- - x : np.ndarray - X values for which the fit function will be evaluated. - minimizer_parameters : dict[str, float] | None, default=None - Dictionary of parameters which will be used in the fit - function. They must be in a dictionary of {'parameter_name': - parameter_value,...}. By default, None. - **kwargs : - Additional arguments. - - Returns - ------- - np.ndarray - Y values calculated at points x for a set of parameters. - - Raises - ------ - TypeError - If ``minimizer_parameters`` is not a dictionary. - """ - if minimizer_parameters is None: - minimizer_parameters = {} - if not isinstance(minimizer_parameters, dict): - raise TypeError('minimizer_parameters must be a dictionary') - - if self._fit_function is None: - # This will also generate self._cached_pars - self._fit_function = self._generate_fit_function() - - minimizer_parameters = self._prepare_parameters(minimizer_parameters) - - return self._fit_function(x, **minimizer_parameters, **kwargs) - def _get_method_kwargs(self, passed_method: str | None = None) -> dict[str, str]: if passed_method is not None: if passed_method not in self.supported_methods(): @@ -215,112 +160,6 @@ def convert_to_par_object(obj): # todo after constraint changes, add type hint: engine Parameter object. """ - def _prepare_parameters(self, parameters: dict[str, float]) -> dict[str, float]: - """ - Prepare the parameters for the minimizer. - - Parameters - ---------- - parameters : dict[str, float] - Dict of parameters for the minimizer with names as keys. - - Returns - ------- - dict[str, float] - Completed parameter dictionary for the minimizer. - """ - pars = self._cached_pars - - for name, item in pars.items(): - parameter_name = MINIMIZER_PARAMETER_PREFIX + str(name) - if parameter_name not in parameters.keys(): - parameters[parameter_name] = item.value - return parameters - - def _generate_fit_function(self) -> Callable: - """ - Using the user supplied ``fit_function``, wrap it in such a way - we can update ``Parameter`` on iterations. - - Returns - ------- - Callable - A fit function which is compatible with bumps models. - """ - # Original fit function - func = self._original_fit_function - # Get a list of `Parameters` - self._cached_pars = {} - self._cached_pars_vals = {} - for parameter in self._object.get_fit_parameters(): - key = parameter.unique_name - self._cached_pars[key] = parameter - self._cached_pars_vals[key] = (parameter.value, parameter.error) - - # Make a new fit function - def _fit_function(x: np.ndarray, **kwargs) -> np.ndarray: - """ - Wrapped fit function which now has an EasyScience compatible - form. - - Parameters - ---------- - x : np.ndarray - Array of data points to be calculated. - **kwargs : - Key word arguments. - - Returns - ------- - np.ndarray - Points calculated at ``x``. - """ - # Update the `Parameter` values and the callback if needed - # TODO THIS IS NOT THREAD SAFE :-( - - for name, value in kwargs.items(): - par_name = name[1:] - if par_name in self._cached_pars.keys(): - # This will take into account constraints - if self._cached_pars[par_name].value != value: - self._cached_pars[par_name].value = value - - # Since we are calling the parameter fset will be called. - # TODO Pre processing here - return_data = func(x) - # TODO Loading or manipulating data here - return return_data - - _fit_function.__signature__ = self._create_signature(self._cached_pars) - return _fit_function - - @staticmethod - def _create_signature(parameters: Dict[int, Parameter]) -> Signature: - """ - Wrap the function signature. - - This is done as lmfit wants the function to be in the form: f = - (x, a=1, b=2)... Where we need to be generic. Note that this - won't hold for much outside of this scope. - """ - wrapped_parameters = [] - wrapped_parameters.append( - InspectParameter('x', InspectParameter.POSITIONAL_OR_KEYWORD, annotation=_empty) - ) - - for name, parameter in parameters.items(): - default_value = parameter.value - - wrapped_parameters.append( - InspectParameter( - MINIMIZER_PARAMETER_PREFIX + str(name), - InspectParameter.POSITIONAL_OR_KEYWORD, - annotation=_empty, - default=default_value, - ) - ) - return Signature(wrapped_parameters) - @staticmethod def _error_from_jacobian( jacobian: np.ndarray, residuals: np.ndarray, confidence: float = 0.95 diff --git a/src/easyscience/fitting/minimizers/minimizer_bumps.py b/src/easyscience/fitting/minimizers/minimizer_bumps.py index 315d856c..8038c2d4 100644 --- a/src/easyscience/fitting/minimizers/minimizer_bumps.py +++ b/src/easyscience/fitting/minimizers/minimizer_bumps.py @@ -3,19 +3,18 @@ from __future__ import annotations -import copy -import math from typing import TYPE_CHECKING from typing import Any from typing import Callable +from typing import cast import numpy as np from bumps.fitters import FIT_AVAILABLE_IDS from bumps.fitters import FITTERS from bumps.fitters import FitDriver -from bumps.names import Curve from bumps.names import FitProblem from bumps.parameter import Parameter as BumpsParameter +from scipy.optimize import OptimizeResult # causes circular import when Parameter is imported # from easyscience.base_classes import ObjBase @@ -24,17 +23,23 @@ from ..available_minimizers import AvailableMinimizers from .bumps_utils import BumpsProgressMonitor from .bumps_utils import EvalCounter +from .bumps_utils import build_curve_problem +from .bumps_utils import parameter_names +from .bumps_utils import parameter_snapshot +from .bumps_utils import to_bumps_parameter +from .bumps_utils import validate_arrays from .minimizer_base import MINIMIZER_PARAMETER_PREFIX from .minimizer_base import MinimizerBase from .utils import FitError from .utils import FitResults if TYPE_CHECKING: - from bumps.dream.state import MCMCDraw + from bumps.fitters import FitBase -FIT_AVAILABLE_IDS_FILTERED = copy.copy(FIT_AVAILABLE_IDS) -# Considered experimental -FIT_AVAILABLE_IDS_FILTERED.remove('pt') +# 'pt' (parallel tempering) is considered experimental and is not exposed. +# Filtered with a comprehension rather than ``list.remove()`` so that importing +# this module does not raise if a future BUMPS release drops the id. +FIT_AVAILABLE_IDS_FILTERED = [fit_id for fit_id in FIT_AVAILABLE_IDS if fit_id != 'pt'] class Bumps(MinimizerBase): @@ -48,10 +53,10 @@ class Bumps(MinimizerBase): def __init__( self, - obj: object, #: ObjBase, + obj: object, fit_function: Callable, minimizer_enum: AvailableMinimizers | None = None, - ): # todo after constraint changes, add type hint: obj: ObjBase # noqa: E501 + ): """ Initialize the fitting engine. @@ -70,7 +75,8 @@ def __init__( @staticmethod def all_methods() -> list[str]: - return FIT_AVAILABLE_IDS_FILTERED + # Copy so callers cannot mutate the module-level list in place. + return list(FIT_AVAILABLE_IDS_FILTERED) @staticmethod def supported_methods() -> list[str]: @@ -88,7 +94,7 @@ def fit( method: str | None = None, tolerance: float | None = None, max_evaluations: int | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, abort_test: Callable[[], bool] | None = None, minimizer_kwargs: dict | None = None, engine_kwargs: dict | None = None, @@ -106,9 +112,14 @@ def fit( weights : np.ndarray Weights for supplied measured points. model : Callable | None, default=None - Optional Model which is being fitted to. By default, None. + Optional BUMPS ``Curve`` which is being fitted to. When + omitted, one is built from ``fit_function`` and the object's + fit parameters. A supplied ``Curve`` must expose ``pars``, + ``x``, ``y`` and ``dy``, since the results are assembled from + them. By default, None. parameters : list[Parameter] | None, default=None - Optional parameters for the fit. By default, None. + Optional parameters for the fit. Ignored when ``model`` is + supplied. By default, None. method : str | None, default=None Method for minimization. By default, None. tolerance : float | None, default=None @@ -118,17 +129,22 @@ def fit( ``steps`` parameter. If ``None``, the default value defined by the selected BUMPS fitter (``fitclass.settings``) is used. By default, None. - progress_callback : Callable[[dict], bool | None] | None, default=None + progress_callback : Callable[[dict], None] | None, default=None Optional callback for progress updates. The payload field - ``iteration`` carries the BUMPS optimizer step index. By - default, None. + ``iteration`` carries the BUMPS optimizer step index. The + return value is ignored — use ``abort_test`` to stop a + running fit. By default, None. abort_test : Callable[[], bool] | None, default=None Optional callback that returns ``True`` to signal that the fit should be aborted. Called periodically during the BUMPS - optimizer iteration loop. + optimizer iteration loop, and once more after the optimizer + returns in order to distinguish an aborted fit from a + converged one, so it must be side-effect free. An aborted fit + returns ``FitResults(success=False)`` rather than raising. minimizer_kwargs : dict | None, default=None Additional keyword arguments passed to the BUMPS minimizer. - By default, None. + The mapping is copied before use, so it is never mutated. By + default, None. engine_kwargs : dict | None, default=None Additional engine keyword arguments. By default, None. **kwargs : Any @@ -137,74 +153,102 @@ def fit( Returns ------- FitResults - Fit results. + Fit results. ``FitResults.iterations`` is the number of BUMPS + *optimizer steps* consumed (the last reported step index plus + one), which is what ``max_evaluations`` budgets against; it is + not comparable to LMFit's ``nfev`` or DFO-LS' ``nf``. The + objective-call count is reported separately as + ``FitResults.n_evaluations``, which is the cross-backend + consistent figure. Note that BUMPS derives the step index from + whatever the selected fitter reports to its monitors, so the + granularity of a "step" varies between fitters. Raises ------ FitError - If the BUMPS fit fails. + If the BUMPS fit raises. A fit that merely fails to converge + is reported as ``FitResults(success=False)`` instead. ValueError - If the input shapes or weights are invalid. + If the input shapes or weights are invalid, or if + ``progress_callback`` is not callable. """ method_dict = self._get_method_kwargs(method) x, y, weights = np.asarray(x), np.asarray(y), np.asarray(weights) - if y.shape != x.shape: - raise ValueError('x and y must have the same shape.') - - if weights.shape != x.shape: - raise ValueError('Weights must have the same shape as x and y.') + validate_arrays(x, y, weights, check_finite_xy=False) - if not np.isfinite(weights).all(): - raise ValueError('Weights cannot be NaN or infinite.') - - if (weights <= 0).any(): - raise ValueError('Weights must be strictly positive and non-zero.') + if progress_callback is not None and not callable(progress_callback): + raise ValueError('progress_callback must be callable') if engine_kwargs is None: engine_kwargs = {} - if minimizer_kwargs is None: - minimizer_kwargs = {} + # Copy rather than mutate: `ftol`/`xtol`/`steps` are injected below, and a + # caller reusing the same mapping for a second fit would otherwise silently + # inherit the settings resolved for the first one. + minimizer_kwargs = {} if minimizer_kwargs is None else dict(minimizer_kwargs) minimizer_kwargs.update(engine_kwargs) method_str = method_dict.get('method', self._method) fitclass = self._resolve_fitclass(method_str) + # Reset the per-fit evaluation counter. A caller-supplied `model` bypasses + # `build_curve_problem`, which is what installs the counter, so without this + # the results would carry the previous fit's objective-call count. + self._eval_counter = None + # Resolve BUMPS-native defaults so the budget reported back to the caller (and # used by the budget-exhaustion check in `_gen_fit_results`) reflects the values # actually consumed by the fitter, even when the caller passes None. + # + # Only values the caller supplied explicitly are pushed back into + # `minimizer_kwargs`. BUMPS pairs an independent `ftol`/`xtol` default per + # fitter (`newton` combines ftol=1e-6 with xtol=1e-12, `amoeba` ftol=1e-8 with + # xtol=1e-6), so collapsing them onto a single resolved value would silently + # tighten the fitter's own convergence criteria on the default path. fitter_settings = dict(fitclass.settings) - if max_evaluations is None: + + if max_evaluations is not None: + minimizer_kwargs['steps'] = max_evaluations + else: max_evaluations = fitter_settings.get('steps') - if tolerance is None: - ftol = fitter_settings.get('ftol') - xtol = fitter_settings.get('xtol') - tols = [t for t in (ftol, xtol) if t is not None] - tolerance = min(tols) if tols else None if tolerance is not None: minimizer_kwargs['ftol'] = tolerance # tolerance for change in function value minimizer_kwargs['xtol'] = ( tolerance # tolerance for change in parameter value, could be an independent value ) - if max_evaluations is not None: - minimizer_kwargs['steps'] = max_evaluations + else: + # Report the stricter of the two BUMPS defaults; nothing is written back. + tols = [ + t + for t in (fitter_settings.get('ftol'), fitter_settings.get('xtol')) + if t is not None + ] + tolerance = min(tols) if tols else None if model is None: - model_function = self._make_model(parameters=parameters) - model = model_function(x, y, weights) + # The Curve comes back directly from the helper: do NOT read it + # from ``problem.fitness``, which is deprecated in BUMPS and warns. + problem, self._eval_counter, model = build_curve_problem( + self, x, y, weights, parameters=parameters + ) + else: + # A caller-supplied model bypasses `build_curve_problem`, which is also + # what populates the parameter cache that `_p_0`, + # `_set_parameter_fit_result` and `_gen_fit_results` all read. Build the + # wrapped fit function here purely for that side effect, so the cache + # describes the current object rather than being empty or left over from + # an earlier fit. + self._fit_function = self._generate_fit_function() + problem = FitProblem(model) self._cached_model = model self._p_0 = {f'p{key}': self._cached_pars[key].value for key in self._cached_pars.keys()} - problem = FitProblem(model) - monitors = [] if progress_callback is not None: - if not callable(progress_callback): - raise ValueError('progress_callback must be callable') monitors.append( BumpsProgressMonitor(problem, progress_callback, self._build_progress_payload) ) @@ -229,25 +273,55 @@ def fit( # Drive the fit through the local FitDriver instance so the supplied # `monitors` (including the optional progress callback monitor) are # invoked. `bumps.fitters.fit` constructs its own driver. - x, fx = driver.fit() - from scipy.optimize import OptimizeResult + # + # Named `best_x` rather than `x` so the caller's independent-variable + # array stays intact for the rest of the method. + best_x, fx = driver.fit() + + # BUMPS signals a failed optimization by returning `None` in place of a + # parameter vector (e.g. Levenberg-Marquardt landing on non-finite + # values); `FitDriver.fit` skips `problem.setp` in that case. Poll + # `abort_test` once more to tell a user-cancelled run apart from a + # converged one, since BUMPS stops quietly either way. + if best_x is None: + success = False + message = 'BUMPS returned no solution; the fit did not converge' + elif abort_test is not None and abort_test(): + success = False + message = 'Fit aborted before convergence' + else: + success = True + message = 'successful termination' # BUMPS' `MonitorRunner.history.step` is populated by the driver itself # (independently of any user-supplied monitors) and exposes the canonical - # last-step index reached by the fitter, so we use it as `nit`. - history_step = getattr(getattr(driver, 'monitor_runner', None), 'history', None) - nit_value = int(history_step.step[0]) if history_step is not None else None + # last-step index reached by the fitter, so we use it as `nit`. `Trace` + # indexes into an internal list, so an empty trace raises `IndexError` + # rather than returning a default — that happens when the fit is aborted + # before the fitter reports its first step. + history = getattr(getattr(driver, 'monitor_runner', None), 'history', None) + step_trace = getattr(history, 'step', None) + nit_value = int(step_trace[0]) if step_trace is not None and len(step_trace) else None + model_results = OptimizeResult( - x=x, - dx=driver.stderr(), + # `driver.stderr()` derives the errors from the covariance at the + # solution, so it cannot be evaluated without one. + x=best_x, + dx=driver.stderr() if best_x is not None else None, fun=fx, - success=True, - status=0, - message='successful termination', + success=success, + status=0 if success else 1, + message=message, nit=nit_value, ) model_results.state = driver.fitter.state - self._set_parameter_fit_result(model_results, stack_status, problem._parameters) + + if best_x is None: + self._restore_parameter_values() + else: + self._set_parameter_fit_result( + model_results, stack_status, parameter_names(problem) + ) results = self._gen_fit_results( model_results, max_evaluations=max_evaluations, @@ -255,27 +329,49 @@ def fit( ) except Exception as e: self._restore_parameter_values() - raise FitError(e) + raise FitError(e) from e finally: global_object.stack.enabled = stack_status return results @staticmethod - def _resolve_fitclass(method: str): - for fitclass in FITTERS: - if fitclass.id == method: - return fitclass - raise FitError(f'Unknown BUMPS fitting method: {method}') + def _resolve_fitclass(method: str) -> type[FitBase]: + """ + Look up the BUMPS fitter class registered under ``method``. + + Parameters + ---------- + method : str + A BUMPS fitter id, e.g. ``'amoeba'``. + + Returns + ------- + type[FitBase] + The matching BUMPS fitter class. + + Raises + ------ + FitError + If no registered fitter carries that id. + """ + # Built per call rather than cached at import time so that fitters + # registered into `FITTERS` after import are still resolvable. + fitclass = {fitclass.id: fitclass for fitclass in FITTERS}.get(method) + if fitclass is None: + raise FitError(f'Unknown BUMPS fitting method: {method}') + # BUMPS annotates `FITTERS` as `List[FitBase]`, but it holds the fitter + # *classes* — `FitDriver` instantiates them as `self.fitclass(problem)`. + return cast('type[FitBase]', fitclass) def _build_progress_payload( - self, problem, iteration: int, point: np.ndarray, nllf: float + self, problem: FitProblem, iteration: int, point: np.ndarray, nllf: float ) -> dict: # Use the nllf already computed by the fitter to avoid a costly # model re-evaluation, and let BUMPS apply its own chisq scaling. chi2 = float(problem.chisq(nllf=nllf, norm=False)) reduced_chi2 = float(problem.chisq(nllf=nllf, norm=True)) - parameter_values = self._current_parameter_snapshot(problem, point) + parameter_values = parameter_snapshot(problem, point) return { 'iteration': iteration, @@ -286,15 +382,6 @@ def _build_progress_payload( 'finished': False, } - def _current_parameter_snapshot(self, problem, point: np.ndarray) -> dict: - labels = problem.labels() - values = problem.getp() if point is None else point - snapshot = {} - for label, value in zip(labels, values): - dict_name = label[len(MINIMIZER_PARAMETER_PREFIX) :] - snapshot[dict_name] = float(value) - return snapshot - def convert_to_pars_obj(self, par_list: list[Parameter] | None = None) -> list[BumpsParameter]: """ Create a container with the ``Parameters`` converted from the @@ -314,10 +401,9 @@ def convert_to_pars_obj(self, par_list: list[Parameter] | None = None) -> list[B if par_list is None: # Assume that we have a ObjBase for which we can obtain a list par_list = self._object.get_fit_parameters() - pars_obj = [self.__class__.convert_to_par_object(obj) for obj in par_list] + pars_obj = [self.convert_to_par_object(obj) for obj in par_list] return pars_obj - # For some reason I have to double staticmethod :-/ @staticmethod def convert_to_par_object(obj: Parameter) -> BumpsParameter: """ @@ -334,386 +420,13 @@ def convert_to_par_object(obj: Parameter) -> BumpsParameter: BumpsParameter Bumps Parameter compatible object. """ - - value = obj.value - - return BumpsParameter( - name=MINIMIZER_PARAMETER_PREFIX + obj.unique_name, - value=value, - bounds=[obj.min, obj.max], - fixed=obj.fixed, - ) - - def _make_model(self, parameters: list[BumpsParameter] | None = None) -> Callable: - """ - Generate a bumps model from the supplied ``fit_function`` and - parameters in the base object. Note that this makes a callable - as it needs to be initialized with *x*, *y*, *weights* - - Weights are converted to dy (standard deviation of y). - - Parameters - ---------- - parameters : list[BumpsParameter] | None, default=None - Optional BUMPS parameters to bind into the model. - - Returns - ------- - Callable - Callable to make a bumps Curve model. - """ - fit_func = EvalCounter(self._generate_fit_function()) - self._eval_counter = fit_func - - def _outer(obj): - - def _make_func(x, y, weights): - bumps_pars = {} - if not parameters: - for name, par in obj._cached_pars.items(): - bumps_pars[MINIMIZER_PARAMETER_PREFIX + str(name)] = ( - obj.convert_to_par_object(par) - ) - else: - for par in parameters: - bumps_pars[MINIMIZER_PARAMETER_PREFIX + par.unique_name] = ( - obj.convert_to_par_object(par) - ) - return Curve(fit_func, x, y, dy=1 / weights, **bumps_pars) - - return _make_func - - return _outer(self) - - def mcmc_sample( - self, - x: np.ndarray, - y: np.ndarray, - weights: np.ndarray, - samples: int = 10000, - burn: int = 2000, - thin: int = 10, - population: int | None = None, - resume_state: MCMCDraw | None = None, - sampler_kwargs: dict | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, - abort_test: Callable[[], bool] | None = None, - ) -> dict: - """ - Run Bayesian MCMC sampling using the BUMPS DREAM sampler. - - Builds a BUMPS ``FitProblem`` from the current model and runs - the DREAM sampler. This is the public minimizer-level entry - point for Bayesian sampling; the higher-level - ``MultiFitter.mcmc_sample`` delegates to this method after - flattening multi-dataset arrays. - - Parameters - ---------- - x : np.ndarray - Flattened independent variable array. - y : np.ndarray - Flattened dependent variable array. - weights : np.ndarray - Flattened weight array. - samples : int, default=10000 - Number of raw samples to draw across all chains, before thinning. - A guaranteed minimum, not an exact count: DREAM advances in - blocks of 10 generations (one generation = one draw per chain) - and stops at the first block boundary at or past ``samples``. - burn : int, default=2000 - Burn-in generations to discard. BUMPS counts ``burn`` in - generations while ``samples`` counts raw draws, so ``burn=500`` - discards ``500 * n_chains`` raw samples. - thin : int, default=10 - Thinning interval — only every ``thin``-th generation is stored. - population : int | None, default=None - BUMPS DREAM population count per parameter (number of parallel - chains): BUMPS creates ``ceil(population * n_parameters)`` chains. - resume_state : MCMCDraw | None, default=None - A BUMPS ``MCMCDraw`` state object from a previous - ``mcmc_sample()`` call (e.g. ``PosteriorResults.sampler_state``). - When provided, DREAM **continues** the saved chain instead of - starting cold. The population, parameter count, and parameter - names must match the current model — a ``ValueError`` is raised - otherwise. - - ``samples`` must be the **total** number of raw samples, not an - increment: to extend an existing chain of ``N`` raw samples by - ``M``, pass ``samples=N + M`` (DREAM keeps only the last - ``samples`` draws in its buffer). The `Sampler.extend` helper - computes this for you. - - ``burn`` is forced to 0 on resume: a previously-converged chain is - never re-burned. - - The ``population`` and ``initializer`` parameters - have **no effect** when ``resume_state`` is provided — they - are determined by the saved state. - - Resuming against *different* data is undefined behaviour (the - chain's likelihood changes underneath it). - sampler_kwargs : dict | None, default=None - Additional keyword arguments forwarded to - ``bumps.fitters.fit``. - progress_callback : Callable[[dict], bool | None] | None, default=None - Optional callback for progress updates during sampling. The - payload dict includes ``iteration`` (DREAM generation - number) and ``sampling: True``. - abort_test : Callable[[], bool] | None, default=None - Optional callback that returns ``True`` to signal that - sampling should be aborted. Called periodically during the - DREAM sampling loop. - - Returns - ------- - dict - Dictionary with keys ``'draws'``, ``'param_names'``, - ``'internal_bumps_object'``, and ``'logp'``. - - Raises - ------ - ValueError - If the input shapes or weights are invalid, if - ``progress_callback`` is not callable, or if ``resume_state`` - is incompatible with the current model (parameter count, - names/order, or population mismatch). - FitError - If DREAM sampling was aborted by the user (via - ``abort_test``). - Exception - Re-raised from DREAM fitting if any unexpected error occurs - (parameter values are restored beforehand). - """ - from bumps.fitters import DreamFit - from bumps.names import FitProblem - - x, y, weights = np.asarray(x), np.asarray(y), np.asarray(weights) - - if not isinstance(samples, int) or samples <= 0: - raise ValueError('samples must be a positive integer.') - if not isinstance(burn, int) or burn < 0: - raise ValueError('burn must be a non-negative integer.') - if not isinstance(thin, int) or thin < 1: - raise ValueError('thin must be a positive integer.') - - if y.shape != x.shape: - raise ValueError('x and y must have the same shape.') - - if not np.isfinite(x).all(): - raise ValueError('x cannot contain NaN or infinite values.') - if not np.isfinite(y).all(): - raise ValueError('y cannot contain NaN or infinite values.') - - if weights.shape != x.shape: - raise ValueError('Weights must have the same shape as x and y.') - - if not np.isfinite(weights).all(): - raise ValueError('Weights cannot be NaN or infinite.') - - if (weights <= 0).any(): - raise ValueError('Weights must be strictly positive and non-zero.') - - # Build the BUMPS Curve model using the minimizer's existing machinery - model_func = self._make_model() - curve = model_func(x, y, weights) - problem = FitProblem(curve) - - pop = population - if resume_state is not None: - pop, burn = self._validate_resume_state(problem, resume_state, population, burn) - - # Build DREAM kwargs. Use the resolved ``pop``, not the raw - # ``population`` argument — on resume ``pop`` is the negative - # absolute chain count that reproduces the saved state's - # population, which BUMPS requires to match. - dream_kwargs: dict = {'samples': samples, 'burn': burn, 'thin': thin} - if pop is not None: - dream_kwargs['pop'] = pop - if sampler_kwargs: - dream_kwargs.update(sampler_kwargs) - - # Build monitors (same pattern as classical Bumps.fit()) - monitors = [] - if progress_callback is not None: - if not callable(progress_callback): - raise ValueError('progress_callback must be callable') - # Compute total DREAM steps for progress display (burn + sampling generations). - # BUMPS DREAM default population count is 10 when not specified by the user. - # A negative ``pop`` (resume) is an absolute chain count. - _dream_default_pop = 10 - pop_val = abs(pop) if pop is not None else _dream_default_pop - _total_steps = burn + (samples + pop_val - 1) // pop_val - monitors.append( - BumpsProgressMonitor( - problem, - progress_callback, - lambda problem, iteration, point, nllf: { - **self._build_sample_progress_payload(problem, iteration, point, nllf), - 'total_steps': _total_steps, - }, - ) - ) - - driver = FitDriver( - fitclass=DreamFit, - problem=problem, - monitors=monitors, - abort_test=abort_test if abort_test is not None else (lambda: False), - **dream_kwargs, - ) - driver.clip() - - from easyscience import global_object - - stack_status = global_object.stack.enabled - global_object.stack.enabled = False - - try: - fit_kwargs = {} - if resume_state is not None: - # Defensive copy: BUMPS mutates the state object in-place - # (via MCMCDraw.resize() — see bumps/dream/core.py allocate_state) - # during resume. Without a copy, the caller's original state - # object is silently altered, making it impossible to compare - # pre- and post-resume state (shape mismatch). See - # https://github.com/easyscience/core/pull/257 - fit_kwargs['fit_state'] = copy.deepcopy(resume_state) - x_opt, fx = driver.fit(**fit_kwargs) - result_state = getattr(driver.fitter, 'state', None) - if result_state is None: - raise FitError('Sampling aborted by user') - except Exception: - self._restore_parameter_values() - raise - finally: - global_object.stack.enabled = stack_status - - _draw = result_state.draw() - draws = _draw.points - param_names = [p.name[len(MINIMIZER_PARAMETER_PREFIX) :] for p in problem._parameters] - logp = _draw.logp - - return { - 'draws': draws, - 'param_names': param_names, - 'internal_bumps_object': result_state, - 'logp': logp, - } - - def _validate_resume_state( - self, - problem: FitProblem, - resume_state: MCMCDraw, - population: int | None, - burn: int, - ) -> tuple[int, int]: - """Check that ``resume_state`` is compatible with ``problem`` and - resolve the population and burn values to use when resuming. - - Parameters - ---------- - problem : FitProblem - The freshly built BUMPS ``FitProblem`` for the current model. - resume_state : MCMCDraw - The saved chain state to resume from. - population : int | None - The caller-supplied population scale factor, or ``None``. - burn : int - The caller-supplied burn-in, ignored (with a warning) on resume. - - Returns - ------- - tuple[int, int] - ``(population, burn)`` to pass to DREAM. The population is - returned as a **negative** number, which BUMPS' - ``initpop.generate`` reads as an absolute chain count, exactly - reproducing the saved state's population. ``burn`` is always 0: - a previously converged chain is never re-burned. - - Raises - ------ - ValueError - If ``resume_state`` is incompatible with the current model - (parameter count, names/order, or population mismatch). - """ - from easyscience import global_object - - logger = global_object.log.getLogger('fitting.bumps') - - # Parameter count - n_params = len(problem._parameters) - if n_params != resume_state.Nvar: - raise ValueError( - f'resume_state has {resume_state.Nvar} parameters but the current ' - f'model has {n_params}. The model must have the same ' - f'number of fitted parameters as when the saved chain was created.' - ) - - prefix = MINIMIZER_PARAMETER_PREFIX - fresh_names = [(p.name or '')[len(prefix) :] for p in problem._parameters] - state_labels = list(resume_state.labels) - if state_labels and all(lbl.startswith(prefix) for lbl in state_labels): - state_names = [lbl[len(prefix) :] for lbl in state_labels] - if fresh_names != state_names: - raise ValueError( - f'Parameter names/order mismatch between the current model ' - f'and resume_state.\n' - f' Current model : {fresh_names}\n' - f' resume_state : {state_names}' - ) - else: - logger.warning( - 'resume_state does not carry parameter names (it was most ' - 'likely reloaded from disk, where BUMPS does not preserve ' - 'labels). Parameter-name validation is skipped; the saved ' - 'chain is matched to the current model by parameter order. ' - 'Ensure this is the same model, with parameters in the same ' - 'order, used to create the chain.' - ) - - # Population. BUMPS creates ``ceil(population * n_params)`` chains - # and requires the resumed state's chain count to match. - if population is not None: - expected_npop = math.ceil(population * n_params) - if expected_npop != resume_state.Npop: - raise ValueError( - f'Requested population ({population}) would produce ' - f'{expected_npop} chains but the saved state has ' - f'{resume_state.Npop} chains. The population cannot ' - f'be changed on resume.' - ) - if burn > 0: - logger.warning( - f'burn={burn} ignored on resume: a previously converged ' - f'chain is not re-burned. Forcing burn=0.' - ) - - # A negative ``pop`` is read by ``bumps.initpop.generate`` as an - # absolute chain count, exactly reproducing the saved population - # without having to recover the original scale factor. - return -int(resume_state.Npop), 0 - - def _build_sample_progress_payload( - self, problem, iteration: int, point: np.ndarray, nllf: float - ) -> dict: - """ - Build a progress payload for Bayesian DREAM sampling steps. - - Called by :class:`BumpsProgressMonitor` at each DREAM - generation. The payload includes ``sampling: True`` so - downstream consumers can distinguish sampling progress from - classical fitting progress. - """ - payload = self._build_progress_payload(problem, iteration, point, nllf) - payload['sampling'] = True - return payload + return to_bumps_parameter(obj) def _set_parameter_fit_result( self, fit_result: Any, stack_status: bool, - par_list: list[BumpsParameter], + par_names: list[str], ) -> None: """ Update parameters to their final values and assign a std error @@ -725,24 +438,28 @@ def _set_parameter_fit_result( BUMPS OptimizeResult containing best-fit values and errors. stack_status : bool Whether the undo stack was enabled. - par_list : list[BumpsParameter] - List of BUMPS parameter objects. + par_names : list[str] + Cached-parameter names in BUMPS problem order, already + stripped of ``MINIMIZER_PARAMETER_PREFIX`` — see + :func:`~easyscience.fitting.minimizers.bumps_utils.parameter_names`. """ from easyscience import global_object pars = self._cached_pars x_result = np.asarray(fit_result.x) - stderr = np.asarray(fit_result.dx) + # Some BUMPS fitters cannot produce a covariance and hand back no errors; + # report those parameters as having no uncertainty rather than failing, + # matching what the LMFit minimizer does when `errorbars` is False. + stderr = None if fit_result.dx is None else np.asarray(fit_result.dx) if stack_status: self._restore_parameter_values() global_object.stack.enabled = True global_object.stack.beginMacro('Fitting routine') - for index, name in enumerate([par.name for par in par_list]): - dict_name = name[len(MINIMIZER_PARAMETER_PREFIX) :] - pars[dict_name].value = x_result[index] - pars[dict_name].error = stderr[index] + for index, name in enumerate(par_names): + pars[name].value = x_result[index] + pars[name].error = 0.0 if stderr is None else stderr[index] if stack_status: global_object.stack.endMacro() @@ -774,8 +491,10 @@ def _gen_fit_results( """ results = FitResults() + # `hasattr`, not a truthiness test: every `FitResults` field starts out + # falsy, so testing the current value would silently discard every kwarg. for name, value in kwargs.items(): - if getattr(results, name, False): + if hasattr(results, name): setattr(results, name, value) n_evaluations = None if self._eval_counter is None else self._eval_counter.count # BUMPS exposes `nit` as the last reported optimizer step index rather than the @@ -808,34 +527,35 @@ def _gen_fit_results( results.p = item results.x = self._cached_model.x results.y_obs = self._cached_model.y + # Costs one extra model evaluation beyond those the optimizer consumed, and + # deliberately so: it runs through the uncounted `self._fit_function`, keeping + # `n_evaluations` a faithful count of optimizer-driven objective calls. results.y_calc = self.evaluate(results.x, minimizer_parameters=results.p) results.y_err = self._cached_model.dy results.n_evaluations = n_evaluations results.iterations = n_steps_used - results.message = '' + # A successful fit carries no message; anything else reports why it stopped. + results.message = ( + '' if fit_results.success else (getattr(fit_results, 'message', '') or '') + ) + if stopped_on_budget: + from easyscience import global_object + results.message = ( f'Fit stopped: reached maximum optimizer steps ({max_evaluations}); ' f'objective evaluated {n_evaluations} times' ) - if stopped_on_budget: - from easyscience import global_object - if tolerance is None: - global_object.log.getLogger('fitting.bumps').warning( - f'Fit did not converge within the maximum optimizer steps of {max_evaluations} ' - f'({n_evaluations} objective evaluations). ' - 'Consider increasing the maximum number of evaluations or adjusting the tolerance.' - ) + reason = 'Fit did not converge within' else: - global_object.log.getLogger('fitting.bumps').warning( - f'Fit did not reach the desired tolerance of {tolerance} within the maximum optimizer steps of {max_evaluations} ' - f'({n_evaluations} objective evaluations). ' - 'Consider increasing the maximum number of evaluations or adjusting the tolerance.' - ) + reason = f'Fit did not reach the desired tolerance of {tolerance} within' + global_object.log.getLogger('fitting.bumps').warning( + f'{reason} the maximum optimizer steps of {max_evaluations} ' + f'({n_evaluations} objective evaluations). ' + 'Consider increasing the maximum number of evaluations or adjusting the tolerance.' + ) - # results.residual = results.y_obs - results.y_calc - # results.goodness_of_fit = np.sum(results.residual**2) results.minimizer_engine = self.__class__ results.fit_args = None results.engine_result = fit_results diff --git a/src/easyscience/fitting/minimizers/minimizer_dfo.py b/src/easyscience/fitting/minimizers/minimizer_dfo.py index 319db849..51c57d63 100644 --- a/src/easyscience/fitting/minimizers/minimizer_dfo.py +++ b/src/easyscience/fitting/minimizers/minimizer_dfo.py @@ -83,7 +83,7 @@ def fit( method: str | None = None, tolerance: float | None = None, max_evaluations: int | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, callback: Callable[[DFOCallbackState], None] | None = None, **kwargs, ) -> FitResults: @@ -108,8 +108,9 @@ def fit( Requested optimizer tolerance. By default, None. max_evaluations : int | None, default=None Maximum number of evaluations. By default, None. - progress_callback : Callable[[dict], bool | None] | None, default=None - Optional callback receiving normalized progress payloads. + progress_callback : Callable[[dict], None] | None, default=None + Optional callback receiving normalized progress payloads. Its + return value is ignored. callback : Callable[[DFOCallbackState], None] | None, default=None Optional native DFO callback. **kwargs : @@ -302,7 +303,7 @@ def wrapped_model(pars_values: List[float]) -> np.ndarray: @staticmethod def _make_progress_adapter( - progress_callback: Callable[[dict], bool | None], + progress_callback: Callable[[dict], None], ) -> Callable[['DFOCallbackState'], None]: """ Create a DFO callback that translates DFOCallbackState into the @@ -310,8 +311,8 @@ def _make_progress_adapter( Parameters ---------- - progress_callback : Callable[[dict], bool | None] - Standard progress callback (dict -> bool|None). + progress_callback : Callable[[dict], None] + Standard progress callback (dict -> None). Returns ------- diff --git a/src/easyscience/fitting/minimizers/minimizer_lmfit.py b/src/easyscience/fitting/minimizers/minimizer_lmfit.py index 31fb543f..30664708 100644 --- a/src/easyscience/fitting/minimizers/minimizer_lmfit.py +++ b/src/easyscience/fitting/minimizers/minimizer_lmfit.py @@ -90,7 +90,7 @@ def fit( method: str | None = None, tolerance: float | None = None, max_evaluations: int | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, minimizer_kwargs: dict | None = None, engine_kwargs: dict | None = None, **kwargs, @@ -116,8 +116,9 @@ def fit( Requested optimizer tolerance. By default, None. max_evaluations : int | None, default=None Maximum number of function evaluations. By default, None. - progress_callback : Callable[[dict], bool | None] | None, default=None - Optional callback receiving normalized progress payloads. + progress_callback : Callable[[dict], None] | None, default=None + Optional callback receiving normalized progress payloads. Its + return value is ignored. minimizer_kwargs : dict | None, default=None Additional keyword arguments passed to LMFit's minimizer. By default, None. @@ -194,7 +195,7 @@ def fit( def _create_iter_callback( self, - progress_callback: Callable[[dict], bool | None] | None, + progress_callback: Callable[[dict], None] | None, ) -> Callable | None: def iter_cb(params, iteration: int, residuals: np.ndarray, *args, **kwargs) -> bool: diff --git a/src/easyscience/fitting/multi_fitter.py b/src/easyscience/fitting/multi_fitter.py index e2d88943..568283c6 100644 --- a/src/easyscience/fitting/multi_fitter.py +++ b/src/easyscience/fitting/multi_fitter.py @@ -56,11 +56,20 @@ def _fit_function_wrapper( Callable Wrapped optimizer function. """ - # Extract of a list of callable functions + # Extract of a list of callable functions. + # ``Fitter._fit_function_wrapper`` reads ``self._fit_function``, so it + # is repointed per dataset inside the loop; the original must be + # restored afterwards or every caller (``Fitter.fit`` aside, which + # snapshots it itself, e.g. sampling) is left with the *last* + # dataset's function on the user-visible ``fit_function`` surface. wrapped_fns = [] - for this_x, this_fun in zip(real_x, self._fit_functions): - self._fit_function = this_fun - wrapped_fns.append(Fitter._fit_function_wrapper(self, this_x, flatten=flatten)) + original_fit_function = self._fit_function + try: + for this_x, this_fun in zip(real_x, self._fit_functions): + self._fit_function = this_fun + wrapped_fns.append(Fitter._fit_function_wrapper(self, this_x, flatten=flatten)) + finally: + self._fit_function = original_fit_function def wrapped_fun(x, **kwargs): # Generate an empty Y based on x diff --git a/src/easyscience/fitting/sampler.py b/src/easyscience/fitting/sampler.py index d3915115..b169c2ad 100644 --- a/src/easyscience/fitting/sampler.py +++ b/src/easyscience/fitting/sampler.py @@ -15,7 +15,7 @@ from easyscience import global_object -from .minimizers.minimizer_base import MINIMIZER_PARAMETER_PREFIX +from .engine_base import PARAMETER_PREFIX if TYPE_CHECKING: # avoid import cycles; only needed for type hints from bumps.dream.state import MCMCDraw @@ -218,9 +218,7 @@ def load_chain(path: str | os.PathLike, skip: int = 0) -> tuple[MCMCDraw, list[s # save_state/load_state does not preserve labels, so a reloaded state # typically carries default labels like ['P0', 'P1', ...]. param_names = [ - lbl[len(MINIMIZER_PARAMETER_PREFIX) :] - if lbl.startswith(MINIMIZER_PARAMETER_PREFIX) - else lbl + lbl[len(PARAMETER_PREFIX) :] if lbl.startswith(PARAMETER_PREFIX) else lbl for lbl in state.labels ] @@ -251,19 +249,9 @@ class SamplingResults: logp: np.ndarray state: MCMCDraw - def to_legacy_dict(self) -> dict: - """Return the legacy dict shape produced by the deprecated - ``mcmc_sample()`` APIs.""" - return { - 'draws': self.draws, - 'param_names': self.param_names, - 'internal_bumps_object': self.state, - 'logp': self.logp, - } - class Sampler: - """Bayesian MCMC sampler for one dataset, backed by a Fitter's BUMPS minimizer. + """Bayesian MCMC sampler for one dataset, backed by the BUMPS DREAM engine. One ``Sampler`` instance represents one chain over one ``(x, y, weights)`` dataset. The data is bound at construction; ``sample()`` and ``extend()`` @@ -274,11 +262,12 @@ class Sampler: effect on the sampler, and there are deliberately no setters: to sample different data, create a new ``Sampler``. - Construct directly with a configured ``Fitter`` (or ``MultiFitter``) whose - minimizer has been switched to ``AvailableMinimizers.Bumps``. **Running a - fit first is not required** — the ``Fitter`` supplies the model and the - minimizer, not a fit result, and sampling from the initial parameter values - works fine. + Construct directly with a configured ``Fitter`` (or ``MultiFitter``). + Sampling is independent of the fitter's minimizer — any minimizer (LMFit, + DFO, BUMPS) may stay active; the only requirement is an installed + ``bumps`` package. **Running a fit first is not required** — the + ``Fitter`` supplies the model and fit function, not a fit result, and + sampling from the initial parameter values works fine. It is often worth fitting first anyway. DREAM seeds its whole starting population inside a tiny ball around the parameters' *current* values @@ -286,14 +275,15 @@ class Sampler: chain in the right region and shortens the burn-in needed to reach the typical set. From a poor initial guess, expect to burn for longer. - The sampler is BUMPS/DREAM-specific for now: the BUMPS check in ``_run()`` - is the seam where another backend would plug in. + The sampler is BUMPS/DREAM-specific for now: the ``DreamSampler`` + construction in ``_run()`` is the single line that knows a concrete + backend exists — the seam where a sampler factory would plug in. Parameters ---------- fitter : Fitter - A configured ``Fitter`` (or ``MultiFitter``) whose minimizer has been - switched to ``AvailableMinimizers.Bumps``. + A configured ``Fitter`` (or ``MultiFitter``) supplying the model and + fit function. Its active minimizer is irrelevant to sampling. x : np.ndarray | list[np.ndarray] Independent variable array (or list of arrays for ``MultiFitter``). y : np.ndarray | list[np.ndarray] @@ -478,52 +468,55 @@ def _run( population: int | None, resume_state: MCMCDraw | None, sampler_kwargs: dict | None, - progress_callback: Callable[[dict], bool | None] | None, + progress_callback: Callable[[dict], None] | None, abort_test: Callable[[], bool] | None, ) -> SamplingResults: """Shared sampling engine for ``sample()`` and ``extend()``. Argument validation for ``samples``/``burn``/``thin`` lives in - ``Bumps.mcmc_sample`` (single source of truth). + ``DreamSampler.run``. """ - # Check the minimizer is BUMPS *before* mutating the fitter — a - # non-BUMPS fitter must not be needlessly rebuilt. - minimizer = self._fitter.minimizer - if not (hasattr(minimizer, 'package') and minimizer.package == 'bumps'): + from .available_minimizers import bumps_engine_available + + if not bumps_engine_available: raise RuntimeError( - 'Bayesian sampling requires a BUMPS minimizer. ' - 'Use ``fitter.switch_minimizer(AvailableMinimizers.Bumps)`` first.' + 'Bayesian sampling requires the bumps package. ' + 'Install it with ``pip install bumps``.' ) + from .samplers.sampler_dream import DreamSampler x_fit, x_new, y_new, w_new, dims = self._fitter._precompute_reshaping( self._x, self._y, self._weights, self._vectorized ) + # Required internal bookkeeping write: MultiFitter's + # ``_fit_function_wrapper`` reads ``_dependent_dims`` to reshape + # multi-dataset output. It is the only fitter attribute sampling + # modifies: the user-visible surface (fit_function, minimizer) is + # never mutated. self._fitter._dependent_dims = dims wrapped = self._fitter._fit_function_wrapper(x_new, flatten=True) merged_kwargs = {**self._default_sampler_kwargs, **(sampler_kwargs or {})} - original_fit_func = self._fitter.fit_function - # Assigning fit_function triggers _update_minimizer() and *rebuilds* - # the minimizer object — it must be re-fetched after this assignment. - self._fitter.fit_function = wrapped - try: - minimizer = self._fitter.minimizer - result = minimizer.mcmc_sample( - x=x_fit, - y=y_new, - weights=w_new, - samples=samples, - burn=burn, - thin=thin, - population=population, - resume_state=resume_state, - sampler_kwargs=merged_kwargs or None, - progress_callback=progress_callback, - abort_test=abort_test, - ) - finally: - self._fitter.fit_function = original_fit_func + # A fresh engine per run is deliberate: it is cheap (the parameter + # cache is built lazily), and per-call construction guarantees the + # chain always sees the fitter's *current* fit function and + # parameter set. Chain continuity lives in ``resume_state``, not in + # engine. + engine = DreamSampler(obj=self._fitter.fit_object, fit_function=wrapped) + result = engine.run( + x=x_fit, + y=y_new, + weights=w_new, + samples=samples, + burn=burn, + thin=thin, + population=population, + resume_state=resume_state, + sampler_kwargs=merged_kwargs or None, + progress_callback=progress_callback, + abort_test=abort_test, + ) results = SamplingResults( draws=result['draws'], @@ -542,7 +535,7 @@ def sample( thin: int = 10, population: int | None = None, sampler_kwargs: dict | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, abort_test: Callable[[], bool] | None = None, ) -> SamplingResults: """Run fresh Bayesian MCMC sampling on the bound data. @@ -571,9 +564,10 @@ def sample( sampler_kwargs : dict | None, default=None Additional keyword arguments forwarded to the BUMPS DREAM sampler (merged over the instance defaults). - progress_callback : Callable[[dict], bool | None] | None, default=None + progress_callback : Callable[[dict], None] | None, default=None Optional callback invoked at each DREAM generation. The payload - dict includes ``iteration`` and ``sampling: True``. + dict includes ``iteration`` and ``sampling: True``. Any return + value is ignored. abort_test : Callable[[], bool] | None, default=None Optional callable that returns ``True`` to abort sampling early. @@ -591,7 +585,7 @@ def sample( Exceptions propagate from the sampling engine: ``ValueError`` if ``samples``, ``burn``, or ``thin`` are invalid, and ``RuntimeError`` - if the active minimizer is not a BUMPS instance. + if the ``bumps`` package is not installed. """ if self._state is not None: global_object.log.getLogger('fitting').warning( @@ -615,7 +609,7 @@ def extend( thin: int = 10, total_samples: int | None = None, sampler_kwargs: dict | None = None, - progress_callback: Callable[[dict], bool | None] | None = None, + progress_callback: Callable[[dict], None] | None = None, abort_test: Callable[[], bool] | None = None, ) -> SamplingResults: """Continue the existing chain with additional samples. @@ -646,8 +640,9 @@ def extend( sampler_kwargs : dict | None, default=None Additional keyword arguments forwarded to the BUMPS DREAM sampler (merged over the instance defaults). - progress_callback : Callable[[dict], bool | None] | None, default=None - Optional callback invoked at each DREAM generation. + progress_callback : Callable[[dict], None] | None, default=None + Optional callback invoked at each DREAM generation. Any return + value is ignored. abort_test : Callable[[], bool] | None, default=None Optional callable that returns ``True`` to abort sampling early. @@ -660,7 +655,8 @@ def extend( ------ RuntimeError If there is no chain to extend (call ``sample()`` or - ``load_state()`` first), or the minimizer is not BUMPS. + ``load_state()`` first), or the ``bumps`` package is not + installed. Notes ----- diff --git a/src/easyscience/fitting/samplers/__init__.py b/src/easyscience/fitting/samplers/__init__.py new file mode 100644 index 00000000..cfdae878 --- /dev/null +++ b/src/easyscience/fitting/samplers/__init__.py @@ -0,0 +1,6 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause + +from .sampler_dream import DreamSampler + +__all__ = ['DreamSampler'] diff --git a/src/easyscience/fitting/samplers/sampler_dream.py b/src/easyscience/fitting/samplers/sampler_dream.py new file mode 100644 index 00000000..b985ad9f --- /dev/null +++ b/src/easyscience/fitting/samplers/sampler_dream.py @@ -0,0 +1,381 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause +"""The BUMPS DREAM MCMC engine — ``DreamSampler``. + +One file per sampling backend, mirroring the one-file-per-minimizer +layout under ``fitting/minimizers/``. When a second MCMC backend +arrives, its ``run()`` signature is formalized as a ``SamplerBase`` ABC +and dispatched via a factory (see discussion easyscience/core#280). +""" + +from __future__ import annotations + +import copy +import math +from typing import TYPE_CHECKING +from typing import Callable + +import numpy as np +from bumps.fitters import FitDriver + +from ..engine_base import PARAMETER_PREFIX +from ..engine_base import EngineBase +from ..minimizers.bumps_utils import BumpsProgressMonitor +from ..minimizers.bumps_utils import build_curve_problem +from ..minimizers.bumps_utils import parameter_names +from ..minimizers.bumps_utils import parameter_snapshot +from ..minimizers.bumps_utils import validate_arrays +from ..minimizers.bumps_utils import validate_run_settings +from ..minimizers.utils import FitError + +if TYPE_CHECKING: + from bumps.dream.state import MCMCDraw + from bumps.names import FitProblem + + +class DreamSampler(EngineBase): + """ + BUMPS DREAM MCMC engine. Runs and resumes chains for one + ``(obj, fit_function)`` binding. + + This is the minimizer-independent home of Bayesian sampling: it + builds its own BUMPS ``FitProblem`` via the shared ``bumps_utils`` + helpers, so sampling no longer requires the ``Fitter``'s active + minimizer to be BUMPS — only an installed ``bumps`` package. + + ``DreamSampler`` is internal machinery; the public entry point is + :class:`easyscience.fitting.Sampler`. + """ + + package = 'bumps' + + def __init__( + self, + obj: object, #: ObjBase, + fit_function: Callable, + ): # todo after constraint changes, add type hint: obj: ObjBase # noqa: E501 + """ + Initialize the sampling engine. + + Parameters + ---------- + obj : object + Object containing the ``Parameter`` instances to sample. + fit_function : Callable + Callable returning model y values for the supplied x values. + """ + super().__init__(obj=obj, fit_function=fit_function) + + def run( + self, + x: np.ndarray, + y: np.ndarray, + weights: np.ndarray | None, + samples: int = 10000, + burn: int = 2000, + thin: int = 10, + population: int | None = None, + resume_state: MCMCDraw | None = None, + sampler_kwargs: dict | None = None, + progress_callback: Callable[[dict], None] | None = None, + abort_test: Callable[[], bool] | None = None, + ) -> dict: + """ + Run Bayesian MCMC sampling using the BUMPS DREAM sampler. + + Builds a BUMPS ``FitProblem`` from the bound object and fit + function and runs the DREAM sampler. This is the engine-level + entry point for Bayesian sampling; the public + :class:`~easyscience.fitting.Sampler` delegates to this method + after flattening its bound data. + + Parameters + ---------- + x : np.ndarray + Flattened independent variable array. + y : np.ndarray + Flattened dependent variable array. + weights : np.ndarray | None + Flattened weight array. Must not be ``None`` — sampling has + no default weighting; a clear ``ValueError`` is raised. + samples : int, default=10000 + Number of raw samples to draw across all chains, before thinning. + A guaranteed minimum, not an exact count: DREAM advances in + blocks of 10 generations (one generation = one draw per chain) + and stops at the first block boundary at or past ``samples``. + burn : int, default=2000 + Burn-in generations to discard. BUMPS counts ``burn`` in + generations while ``samples`` counts raw draws, so ``burn=500`` + discards ``500 * n_chains`` raw samples. + thin : int, default=10 + Thinning interval — only every ``thin``-th generation is stored. + population : int | None, default=None + BUMPS DREAM population count per parameter (number of parallel + chains): BUMPS creates ``ceil(population * n_parameters)`` chains. + resume_state : MCMCDraw | None, default=None + A BUMPS ``MCMCDraw`` state object from a previous ``run()`` + call. When provided, DREAM **continues** the saved chain + instead of starting cold. The population, parameter count, + and parameter names must match the current model — a + ``ValueError`` is raised otherwise. + + ``samples`` must be the **total** number of raw samples, not an + increment: to extend an existing chain of ``N`` raw samples by + ``M``, pass ``samples=N + M`` (DREAM keeps only the last + ``samples`` draws in its buffer). The `Sampler.extend` helper + computes this for you. + + ``burn`` is forced to 0 on resume: a previously-converged chain is + never re-burned. + + The ``population`` and ``initializer`` parameters + have **no effect** when ``resume_state`` is provided — they + are determined by the saved state. + + Resuming against *different* data is undefined behaviour (the + chain's likelihood changes underneath it). + sampler_kwargs : dict | None, default=None + Additional keyword arguments forwarded to + ``bumps.fitters.fit``. + progress_callback : Callable[[dict], None] | None, default=None + Optional callback for progress updates during sampling. The + payload dict includes ``iteration`` (DREAM generation + number) and ``sampling: True``. Any return value is ignored. + abort_test : Callable[[], bool] | None, default=None + Optional callback that returns ``True`` to signal that + sampling should be aborted. Called periodically during the + DREAM sampling loop. + + Returns + ------- + dict + Dictionary with keys ``'draws'``, ``'param_names'``, + ``'internal_bumps_object'``, and ``'logp'``. + + Raises + ------ + ValueError + If the input shapes or weights are invalid, if + ``progress_callback`` is not callable, or if ``resume_state`` + is incompatible with the current model (parameter count, + names/order, or population mismatch). + FitError + If DREAM sampling was aborted by the user (via + ``abort_test``). + Exception + Re-raised from DREAM fitting if any unexpected error occurs + (parameter values are restored beforehand). + """ + from bumps.fitters import DreamFit + + if weights is None: + raise ValueError( + 'weights must not be None for Bayesian sampling. Pass ' + 'measurement weights (e.g. ``1 / sigma``) matching x and y.' + ) + x, y, weights = np.asarray(x), np.asarray(y), np.asarray(weights) + + validate_run_settings(samples, burn, thin) + validate_arrays(x, y, weights, check_finite_xy=True) + + # Build the BUMPS Curve model around the engine's wrapped fit function + problem, _, _ = build_curve_problem(self, x, y, weights) + + pop = population + if resume_state is not None: + pop, burn = self._validate_resume_state(problem, resume_state, population, burn) + + # Build DREAM kwargs. Use the resolved ``pop``, not the raw + # ``population`` argument — on resume ``pop`` is the negative + # absolute chain count that reproduces the saved state's + # population, which BUMPS requires to match. + dream_kwargs: dict = {'samples': samples, 'burn': burn, 'thin': thin} + if pop is not None: + dream_kwargs['pop'] = pop + if sampler_kwargs: + dream_kwargs.update(sampler_kwargs) + + # Build monitors (same pattern as classical Bumps.fit()) + monitors = [] + if progress_callback is not None: + if not callable(progress_callback): + raise ValueError('progress_callback must be callable') + # Compute total DREAM steps for progress display (burn + sampling generations). + # BUMPS DREAM default population count is 10 when not specified by the user. + # A negative ``pop`` (resume) is an absolute chain count. + _dream_default_pop = 10 + pop_val = abs(pop) if pop is not None else _dream_default_pop + _total_steps = burn + (samples + pop_val - 1) // pop_val + monitors.append( + BumpsProgressMonitor( + problem, + progress_callback, + lambda problem, iteration, point, nllf: { + **self._build_sample_progress_payload(problem, iteration, point, nllf), + 'total_steps': _total_steps, + }, + ) + ) + + driver = FitDriver( + fitclass=DreamFit, + problem=problem, + monitors=monitors, + abort_test=abort_test if abort_test is not None else (lambda: False), + **dream_kwargs, + ) + driver.clip() + + from easyscience import global_object + + stack_status = global_object.stack.enabled + global_object.stack.enabled = False + + try: + fit_kwargs = {} + if resume_state is not None: + # Defensive copy: BUMPS mutates the state object in-place + # (via MCMCDraw.resize() — see bumps/dream/core.py allocate_state) + # during resume. Without a copy, the caller's original state + # object is silently altered, making it impossible to compare + # pre- and post-resume state (shape mismatch). See + # https://github.com/easyscience/core/pull/257 + fit_kwargs['fit_state'] = copy.deepcopy(resume_state) + x_opt, fx = driver.fit(**fit_kwargs) + result_state = getattr(driver.fitter, 'state', None) + if result_state is None: + raise FitError('Sampling aborted by user') + except Exception: + self._restore_parameter_values() + raise + finally: + global_object.stack.enabled = stack_status + + _draw = result_state.draw() + + return { + 'draws': _draw.points, + 'param_names': parameter_names(problem), + 'internal_bumps_object': result_state, + 'logp': _draw.logp, + } + + def _validate_resume_state( + self, + problem: FitProblem, + resume_state: MCMCDraw, + population: int | None, + burn: int, + ) -> tuple[int, int]: + """Check that ``resume_state`` is compatible with ``problem`` and + resolve the population and burn values to use when resuming. + + Parameters + ---------- + problem : FitProblem + The freshly built BUMPS ``FitProblem`` for the current model. + resume_state : MCMCDraw + The saved chain state to resume from. + population : int | None + The caller-supplied population scale factor, or ``None``. + burn : int + The caller-supplied burn-in, ignored (with a warning) on resume. + + Returns + ------- + tuple[int, int] + ``(population, burn)`` to pass to DREAM. The population is + returned as a **negative** number, which BUMPS' + ``initpop.generate`` reads as an absolute chain count, exactly + reproducing the saved state's population. ``burn`` is always 0: + a previously converged chain is never re-burned. + + Raises + ------ + ValueError + If ``resume_state`` is incompatible with the current model + (parameter count, names/order, or population mismatch). + """ + from easyscience import global_object + + logger = global_object.log.getLogger('fitting.bumps') + + # Parameter count + n_params = len(problem._parameters) + if n_params != resume_state.Nvar: + raise ValueError( + f'resume_state has {resume_state.Nvar} parameters but the current ' + f'model has {n_params}. The model must have the same ' + f'number of fitted parameters as when the saved chain was created.' + ) + + prefix = PARAMETER_PREFIX + fresh_names = [(p.name or '')[len(prefix) :] for p in problem._parameters] + state_labels = list(resume_state.labels) + if state_labels and all(lbl.startswith(prefix) for lbl in state_labels): + state_names = [lbl[len(prefix) :] for lbl in state_labels] + if fresh_names != state_names: + raise ValueError( + f'Parameter names/order mismatch between the current model ' + f'and resume_state.\n' + f' Current model : {fresh_names}\n' + f' resume_state : {state_names}' + ) + else: + logger.warning( + 'resume_state does not carry parameter names (it was most ' + 'likely reloaded from disk, where BUMPS does not preserve ' + 'labels). Parameter-name validation is skipped; the saved ' + 'chain is matched to the current model by parameter order. ' + 'Ensure this is the same model, with parameters in the same ' + 'order, used to create the chain.' + ) + + # Population. BUMPS creates ``ceil(population * n_params)`` chains + # and requires the resumed state's chain count to match. + if population is not None: + expected_npop = math.ceil(population * n_params) + if expected_npop != resume_state.Npop: + raise ValueError( + f'Requested population ({population}) would produce ' + f'{expected_npop} chains but the saved state has ' + f'{resume_state.Npop} chains. The population cannot ' + f'be changed on resume.' + ) + if burn > 0: + logger.warning( + f'burn={burn} ignored on resume: a previously converged ' + f'chain is not re-burned. Forcing burn=0.' + ) + + # A negative ``pop`` is read by ``bumps.initpop.generate`` as an + # absolute chain count, exactly reproducing the saved population + # without having to recover the original scale factor. + return -int(resume_state.Npop), 0 + + def _build_sample_progress_payload( + self, problem, iteration: int, point: np.ndarray, nllf: float + ) -> dict: + """ + Build a progress payload for Bayesian DREAM sampling steps. + + Called by :class:`BumpsProgressMonitor` at each DREAM + generation. The payload includes ``sampling: True`` so + downstream consumers can distinguish sampling progress from + classical fitting progress; the remaining keys match the + classical-fit payload built by the minimizers. + """ + # Use the nllf already computed by the sampler to avoid a costly + # model re-evaluation, and let BUMPS apply its own chisq scaling. + chi2 = float(problem.chisq(nllf=nllf, norm=False)) + reduced_chi2 = float(problem.chisq(nllf=nllf, norm=True)) + + return { + 'iteration': iteration, + 'chi2': chi2, + 'reduced_chi2': reduced_chi2, + 'parameter_values': parameter_snapshot(problem, point), + 'refresh_plots': False, + 'finished': False, + 'sampling': True, + } diff --git a/tests/integration/fitting/test_fitter.py b/tests/integration/fitting/test_fitter.py index df17a65c..5224da73 100644 --- a/tests/integration/fitting/test_fitter.py +++ b/tests/integration/fitting/test_fitter.py @@ -354,6 +354,35 @@ def test_bumps_methods(fit_method): check_fit_results(result, sp_sin, ref_sin, x) +@pytest.mark.fast +def test_bumps_fit_emits_no_fitness_deprecation_warning(): + """Regression (CR-1): the classical BUMPS fit path must not read the + deprecated ``FitProblem.fitness`` property, which emits a ``UserWarning`` + on bumps >= 1.0.4 — the ``Curve`` comes back from ``build_curve_problem`` + directly.""" + import warnings + + ref_sin = AbsSin(0.2, np.pi) + sp_sin = AbsSin(0.354, 3.05) + + x = np.linspace(0, 5, 200) + weights = np.ones_like(x) + y = ref_sin(x) + + sp_sin.offset.fixed = False + sp_sin.phase.fixed = False + + f = Fitter(sp_sin, sp_sin) + f.switch_minimizer('Bumps') + + with warnings.catch_warnings(record=True) as caught: + warnings.simplefilter('always') + f.fit(x, y, weights=weights) + + fitness_warnings = [w for w in caught if 'fitness' in str(w.message)] + assert fitness_warnings == [] + + @pytest.mark.fast @pytest.mark.parametrize( 'fit_engine', diff --git a/tests/integration/fitting/test_sampler.py b/tests/integration/fitting/test_sampler.py index 032cc39f..a0a0a217 100644 --- a/tests/integration/fitting/test_sampler.py +++ b/tests/integration/fitting/test_sampler.py @@ -62,8 +62,13 @@ def __call__(self, x): ) -def _bumps_fitter_and_data(): - """Build a 2-parameter BUMPS MultiFitter over a small sine model.""" +def _fitter_and_data(): + """Build a 2-parameter MultiFitter over a small sine model. + + The fitter keeps its default (LMFit) minimizer: sampling no longer + requires switching to BUMPS, only an installed ``bumps`` package. + """ + pytest.importorskip('bumps') ref_sin = AbsSin(0.2, np.pi) sp = AbsSin(0.354, 3.05) sp.offset.fixed = False @@ -72,10 +77,6 @@ def _bumps_fitter_and_data(): y = ref_sin(x) weights = np.ones_like(x) f = MultiFitter([sp], [sp]) - try: - f.switch_minimizer('Bumps') - except AttributeError: - pytest.skip('BUMPS is not installed') return f, sp, x, y, weights @@ -84,8 +85,8 @@ class TestSampler: @pytest.mark.filterwarnings('ignore::UserWarning') def test_sample_returns_results_object(self): - """sample() returns a populated SamplingResults; to_legacy_dict() has the legacy shape.""" - f, sp, x, y, weights = _bumps_fitter_and_data() + """sample() returns a populated SamplingResults, cached on the sampler.""" + f, sp, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) results = sampler.sample(samples=100, burn=20, thin=2) @@ -106,11 +107,6 @@ def test_sample_returns_results_object(self): assert sampler.draws is results.draws assert sampler.param_names == results.param_names - # legacy dict shape - legacy = results.to_legacy_dict() - assert set(legacy.keys()) == {'draws', 'param_names', 'internal_bumps_object', 'logp'} - assert legacy['internal_bumps_object'] is results.state - @pytest.mark.filterwarnings('ignore::UserWarning') def test_sample_multi_dataset(self): """Multi-dataset sampling via Sampler(f, ...) has correct param_names.""" @@ -133,11 +129,8 @@ def test_sample_multi_dataset(self): sp_sin_1.phase.fixed = False sp_line.c.fixed = False + pytest.importorskip('bumps') f = MultiFitter([sp_sin_1, sp_line], [sp_sin_1, sp_line]) - try: - f.switch_minimizer('Bumps') - except AttributeError: - pytest.skip('BUMPS is not installed') sampler = Sampler(f, [x1, x2], [y1, y2], [weights, weights]) results = sampler.sample(samples=100, burn=20, thin=2) @@ -149,7 +142,7 @@ def test_sample_multi_dataset(self): def test_sample_population(self): """Passing population should succeed and produce valid draws.""" - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) results = sampler.sample(samples=100, burn=20, thin=2, population=5) @@ -169,11 +162,8 @@ def test_sample_vectorized_2d(self): sp.offset.fixed = False sp.phase.fixed = False + pytest.importorskip('bumps') f = MultiFitter([sp], [sp]) - try: - f.switch_minimizer('Bumps') - except AttributeError: - pytest.skip('BUMPS is not installed') sampler = Sampler(f, [x2D], [y2D], [weights], vectorized=True) results = sampler.sample(samples=100, burn=20, thin=2) @@ -185,17 +175,45 @@ def test_sample_vectorized_2d(self): @pytest.mark.filterwarnings('ignore::UserWarning') def test_fit_function_restored_on_success(self): """fit_function must be restored after a successful sample().""" - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) original_func = f.fit_function sampler.sample(samples=100, burn=20, thin=2) assert f.fit_function is original_func + @pytest.mark.filterwarnings('ignore::UserWarning') + def test_fit_function_untouched_multi_dataset(self): + """With 2+ datasets the per-dataset wrapping in MultiFitter must not + leave fit_function pointing at the LAST dataset's function after + sampling (regression: the single-dataset variant above is vacuous for + this bug because last == first == original).""" + ref_sin = AbsSin(0.2, np.pi) + sp_sin = AbsSin(0.354, 3.05) + sp_line = Line(0.43, 6.1) + sp_sin.offset.fixed = False + sp_line.c.fixed = False + + x1 = np.linspace(0, 5, 50) + y1 = ref_sin(x1) + x2 = np.copy(x1) + y2 = Line(1, 4.6)(x2) + weights = np.ones_like(x1) + + pytest.importorskip('bumps') + f = MultiFitter([sp_sin, sp_line], [sp_sin, sp_line]) + original = f.fit_function + assert original is sp_sin # two distinct per-dataset functions + + sampler = Sampler(f, [x1, x2], [y1, y2], [weights, weights]) + sampler.sample(samples=50, burn=5, thin=1) + + assert f.fit_function is original + @pytest.mark.filterwarnings('ignore::UserWarning') def test_sampler_kwargs_forwarded(self): """Per-call sampler_kwargs dict is forwarded to the BUMPS DREAM sampler.""" - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) results = sampler.sample(samples=100, burn=20, thin=2, sampler_kwargs={'init': 'random'}) @@ -204,33 +222,47 @@ def test_sampler_kwargs_forwarded(self): assert results.draws.shape[0] > 0 @pytest.mark.filterwarnings('ignore::UserWarning') - def test_default_sampler_kwargs_merged(self): + def test_default_sampler_kwargs_merged(self, monkeypatch): """Constructor-level sampler_kwargs defaults are used; per-call kwargs win.""" - f, _, x, y, weights = _bumps_fitter_and_data() + from easyscience.fitting.samplers.sampler_dream import DreamSampler + + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights], sampler_kwargs={'init': 'random'}) captured = {} - original_mcmc_sample = type(f.minimizer).mcmc_sample + original_run = DreamSampler.run def spy(self, **kwargs): captured.update(kwargs.get('sampler_kwargs') or {}) - return original_mcmc_sample(self, **kwargs) + return original_run(self, **kwargs) + + monkeypatch.setattr(DreamSampler, 'run', spy) + + sampler.sample(samples=100, burn=20, thin=2) + assert captured == {'init': 'random'} - try: - type(f.minimizer).mcmc_sample = spy - sampler.sample(samples=100, burn=20, thin=2) - assert captured == {'init': 'random'} + captured.clear() + sampler.sample(samples=100, burn=20, thin=2, sampler_kwargs={'init': 'lhs'}) + assert captured == {'init': 'lhs'} # per-call overrides default - captured.clear() - sampler.sample(samples=100, burn=20, thin=2, sampler_kwargs={'init': 'lhs'}) - assert captured == {'init': 'lhs'} # per-call overrides default - finally: - type(f.minimizer).mcmc_sample = original_mcmc_sample + @pytest.mark.filterwarnings('ignore::UserWarning') + def test_sample_with_lmfit_minimizer_active(self): + """Sampling works without switching the fitter's minimizer to BUMPS — + the new capability enabled by the ``DreamSampler`` engine (#280).""" + f, _, x, y, weights = _fitter_and_data() + assert f.minimizer.package == 'lmfit' # the default LMFit minimizer + + sampler = Sampler(f, [x], [y], [weights]) + results = sampler.sample(samples=100, burn=20, thin=2) + + assert results.draws.shape[0] > 0 + # The active minimizer is untouched by sampling. + assert f.minimizer.package == 'lmfit' @pytest.mark.filterwarnings('ignore::UserWarning') def test_extend_chain(self): """extend(additional_samples=) continues the chain; ring-buffer math is done for the user.""" - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) first = sampler.sample(samples=100, burn=20, thin=1) @@ -250,7 +282,7 @@ def test_extend_with_thinning_keeps_existing_draws(self): generations (``Ngen * Npop``), not from the retained-draw count, which BUMPS divides by the thinning interval. """ - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) first = sampler.sample(samples=1000, burn=20, thin=10) @@ -264,7 +296,7 @@ def test_extend_with_thinning_keeps_existing_draws(self): @pytest.mark.filterwarnings('ignore::UserWarning') def test_extend_total_samples_override(self): """extend(total_samples=) bypasses the additional_samples arithmetic.""" - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) sampler.sample(samples=100, burn=20, thin=1) @@ -282,7 +314,7 @@ def test_extend_after_save_load_roundtrip(self, tmp_path, caplog): """ import logging - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) first = sampler.sample(samples=100, burn=20, thin=1) @@ -308,7 +340,7 @@ def test_extend_preserves_nondefault_population(self): saved state on resume, otherwise BUMPS regenerates the default population and raises ``Cannot change Nvar, Npop or Ncr on resize``. """ - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) first = sampler.sample(samples=100, burn=20, thin=1, population=5) @@ -327,7 +359,7 @@ def test_save_warns_when_fingerprint_unavailable(self, tmp_path, caplog, monkeyp logs a warning and records ``null`` in the sidecar.""" import logging - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) sampler.sample(samples=100, burn=20, thin=2) @@ -344,7 +376,7 @@ def test_save_warns_when_fingerprint_unavailable(self, tmp_path, caplog, monkeyp @pytest.mark.filterwarnings('ignore::UserWarning') def test_load_state_populates_results(self, tmp_path): """A freshly loaded sampler reports draws/logp/param_names without resampling.""" - f, sp, x, y, weights = _bumps_fitter_and_data() + f, sp, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) first = sampler.sample(samples=100, burn=20, thin=2) @@ -377,7 +409,7 @@ def test_load_short_chain_regression(self, tmp_path): reader collapses it to a 1-D array and ``load_state`` raises ``IndexError`` without the 2-D coercion workaround. """ - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) sampler.sample(samples=20, burn=5, thin=1) @@ -391,7 +423,7 @@ def test_load_short_chain_regression(self, tmp_path): @pytest.mark.filterwarnings('ignore::UserWarning') def test_load_fingerprint_mismatch_warns(self, tmp_path, caplog): """Loading a chain into a sampler bound to different data warns.""" - f, _, x, y, weights = _bumps_fitter_and_data() + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) sampler.sample(samples=100, burn=20, thin=2) diff --git a/tests/unit/fitting/minimizers/bumps_utils/__init__.py b/tests/unit/fitting/minimizers/bumps_utils/__init__.py new file mode 100644 index 00000000..4e798e20 --- /dev/null +++ b/tests/unit/fitting/minimizers/bumps_utils/__init__.py @@ -0,0 +1,2 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause diff --git a/tests/unit/fitting/minimizers/bumps_utils/test_problem.py b/tests/unit/fitting/minimizers/bumps_utils/test_problem.py new file mode 100644 index 00000000..affe454d --- /dev/null +++ b/tests/unit/fitting/minimizers/bumps_utils/test_problem.py @@ -0,0 +1,179 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause +"""Unit tests for the shared BUMPS problem-construction helpers.""" + +from unittest.mock import MagicMock + +import numpy as np +import pytest + +import easyscience.fitting.minimizers.bumps_utils.problem +from easyscience.fitting.minimizers.bumps_utils import build_curve_problem +from easyscience.fitting.minimizers.bumps_utils import parameter_names +from easyscience.fitting.minimizers.bumps_utils import parameter_snapshot +from easyscience.fitting.minimizers.bumps_utils import to_bumps_parameter + + +class TestToBumpsParameter: + def test_convert_parameter_object(self) -> None: + from easyscience.variable import Parameter + + param = Parameter('thickness', 42.0, min=0.0, max=100.0) + param.fixed = False + + result = to_bumps_parameter(param) + + # to_bumps_parameter uses obj.unique_name which is auto-assigned + assert result.name.startswith('p') + assert result.value == 42.0 + assert result.bounds == (0.0, 100.0) + assert result.fixed is False + + def test_convert_fixed_parameter(self) -> None: + from easyscience.variable import Parameter + + param = Parameter('roughness', 5.0, min=0.0, max=20.0) + param.fixed = True + + result = to_bumps_parameter(param) + + assert result.name.startswith('p') + assert result.fixed is True + + +class TestBuildCurveProblem: + """Curve/FitProblem assembly, with the BUMPS classes mocked out.""" + + @pytest.fixture(autouse=True) + def _mock_bumps_classes(self, monkeypatch): + self.mock_curve_cls = MagicMock(return_value='curve') + self.mock_problem_cls = MagicMock(return_value='problem') + monkeypatch.setattr( + easyscience.fitting.minimizers.bumps_utils.problem, 'Curve', self.mock_curve_cls + ) + monkeypatch.setattr( + easyscience.fitting.minimizers.bumps_utils.problem, + 'FitProblem', + self.mock_problem_cls, + ) + self.mock_convert = MagicMock(side_effect=lambda par: f'converted-{par.unique_name}') + monkeypatch.setattr( + easyscience.fitting.minimizers.bumps_utils.problem, + 'to_bumps_parameter', + self.mock_convert, + ) + + @staticmethod + def _engine_with_cached_pars(cached_pars): + engine = MagicMock() + engine._generate_fit_function = MagicMock( + return_value=MagicMock(return_value=np.array([2.0])) + ) + engine._cached_pars = cached_pars + return engine + + def test_uses_cached_parameters_by_default(self): + cached_par = MagicMock() + cached_par.unique_name = 'alpha' + engine = self._engine_with_cached_pars({'alpha': cached_par}) + + problem, counter, curve = build_curve_problem( + engine, np.array([1.0]), np.array([2.0]), np.array([4.0]) + ) + + assert problem == 'problem' + # The Curve is surfaced directly so callers never have to read the + # deprecated ``FitProblem.fitness`` property (CR-1). + assert curve == 'curve' + engine._generate_fit_function.assert_called_once_with() + self.mock_convert.assert_called_once_with(cached_par) + assert self.mock_curve_cls.call_args.kwargs['palpha'] == 'converted-alpha' + self.mock_problem_cls.assert_called_once_with('curve') + + def test_explicit_parameters_override_cache(self): + engine = self._engine_with_cached_pars({'alpha': MagicMock(unique_name='alpha')}) + explicit = MagicMock() + explicit.unique_name = 'beta' + + build_curve_problem( + engine, np.array([1.0]), np.array([2.0]), np.array([4.0]), parameters=[explicit] + ) + + self.mock_convert.assert_called_once_with(explicit) + assert 'pbeta' in self.mock_curve_cls.call_args.kwargs + assert 'palpha' not in self.mock_curve_cls.call_args.kwargs + + def test_curve_receives_data_and_dy(self): + """weights are converted to dy = 1 / weights.""" + engine = self._engine_with_cached_pars({}) + x = np.array([1.0, 2.0]) + y = np.array([10.0, 20.0]) + weights = np.array([2.0, 4.0]) + + build_curve_problem(engine, x, y, weights) + + call = self.mock_curve_cls.call_args + np.testing.assert_array_equal(call.args[1], x) + np.testing.assert_array_equal(call.args[2], y) + np.testing.assert_array_equal(call.kwargs['dy'], 1 / weights) + + def test_counter_wraps_fit_function(self): + """The returned EvalCounter wraps the wrapped fit function and counts calls.""" + inner = MagicMock(return_value=np.array([11.0, 22.0])) + engine = self._engine_with_cached_pars({}) + engine._generate_fit_function = MagicMock(return_value=inner) + + _, counter, _ = build_curve_problem( + engine, np.array([1.0]), np.array([2.0]), np.array([4.0]) + ) + + # The counter itself is what Curve receives as the fit function. + assert self.mock_curve_cls.call_args.args[0] is counter + assert counter.count == 0 + counter(np.array([1.0])) + assert counter.count == 1 + inner.assert_called_once() + + +class TestParameterNames: + def test_strips_prefix(self): + params = [] + for name in ('palpha', 'pbeta'): + p = MagicMock() + p.name = name + params.append(p) + problem = MagicMock() + problem._parameters = params + + assert parameter_names(problem) == ['alpha', 'beta'] + + def test_tolerates_none_name(self): + p = MagicMock() + p.name = None + problem = MagicMock() + problem._parameters = [p] + + assert parameter_names(problem) == [''] + + +class TestParameterSnapshot: + def test_snapshot_from_point(self) -> None: + mock_problem = MagicMock() + mock_problem.labels.return_value = ['palpha', 'pbeta'] + + point = np.array([1.5, 2.5]) + + snapshot = parameter_snapshot(mock_problem, point) + + assert snapshot == {'alpha': 1.5, 'beta': 2.5} + mock_problem.getp.assert_not_called() + + def test_snapshot_falls_back_to_getp(self) -> None: + mock_problem = MagicMock() + mock_problem.labels.return_value = ['palpha'] + mock_problem.getp.return_value = np.array([3.5]) + + snapshot = parameter_snapshot(mock_problem, None) + + assert snapshot == {'alpha': 3.5} + mock_problem.getp.assert_called_once() diff --git a/tests/unit/fitting/minimizers/bumps_utils/test_validation.py b/tests/unit/fitting/minimizers/bumps_utils/test_validation.py new file mode 100644 index 00000000..c69d86f8 --- /dev/null +++ b/tests/unit/fitting/minimizers/bumps_utils/test_validation.py @@ -0,0 +1,88 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause +"""Unit tests for the shared BUMPS input validation helpers.""" + +import numpy as np +import pytest + +from easyscience.fitting.minimizers.bumps_utils import validate_arrays +from easyscience.fitting.minimizers.bumps_utils import validate_run_settings + + +class TestValidateRunSettings: + @pytest.mark.parametrize( + 'kwargs, match', + [ + ({'samples': 0}, 'samples must be a positive integer'), + ({'samples': -1}, 'samples must be a positive integer'), + ({'samples': 10.0}, 'samples must be a positive integer'), + # bool is an int subclass; True must not sneak through as 1. + ({'samples': True}, 'samples must be a positive integer'), + ({'burn': -1}, 'burn must be a non-negative integer'), + ({'burn': 1.5}, 'burn must be a non-negative integer'), + ({'burn': True}, 'burn must be a non-negative integer'), + ({'burn': False}, 'burn must be a non-negative integer'), + ({'thin': 0}, 'thin must be a positive integer'), + ({'thin': 2.0}, 'thin must be a positive integer'), + ({'thin': True}, 'thin must be a positive integer'), + ], + ) + def test_invalid_settings_raise(self, kwargs, match): + settings = {'samples': 10, 'burn': 0, 'thin': 1} + settings.update(kwargs) + with pytest.raises(ValueError, match=match): + validate_run_settings(**settings) + + def test_valid_settings_pass(self): + validate_run_settings(samples=1, burn=0, thin=1) + validate_run_settings(samples=10000, burn=2000, thin=10) + + +class TestValidateArrays: + @staticmethod + def _data(): + return { + 'x': np.array([1.0, 2.0]), + 'y': np.array([0.1, 0.2]), + 'weights': np.array([1.0, 1.0]), + } + + @pytest.mark.parametrize( + 'overrides, match', + [ + ({'y': np.array([0.1])}, 'x and y must have the same shape'), + ({'weights': np.array([1.0])}, 'Weights must have the same shape'), + ({'weights': np.array([1.0, np.nan])}, 'Weights cannot be NaN'), + ({'weights': np.array([1.0, np.inf])}, 'Weights cannot be NaN'), + ({'weights': np.array([1.0, 0.0])}, 'Weights must be strictly positive'), + ({'weights': np.array([1.0, -1.0])}, 'Weights must be strictly positive'), + ], + ) + @pytest.mark.parametrize('check_finite_xy', [True, False]) + def test_shared_checks_raise(self, overrides, match, check_finite_xy): + """Shape and weight checks apply on both the fit and sampling paths.""" + data = self._data() + data.update(overrides) + with pytest.raises(ValueError, match=match): + validate_arrays(**data, check_finite_xy=check_finite_xy) + + @pytest.mark.parametrize( + 'overrides, match', + [ + ({'x': np.array([1.0, np.nan])}, 'x cannot contain NaN'), + ({'x': np.array([1.0, np.inf])}, 'x cannot contain NaN'), + ({'y': np.array([0.1, np.nan])}, 'y cannot contain NaN'), + ({'y': np.array([0.1, np.inf])}, 'y cannot contain NaN'), + ], + ) + def test_finite_xy_checked_only_when_requested(self, overrides, match): + """x/y finiteness is enforced for sampling but not for the classical + fit path, preserving the fit path's historically permissive behaviour.""" + data = self._data() + data.update(overrides) + with pytest.raises(ValueError, match=match): + validate_arrays(**data, check_finite_xy=True) + validate_arrays(**data, check_finite_xy=False) # must not raise + + def test_valid_arrays_pass(self): + validate_arrays(**self._data(), check_finite_xy=True) diff --git a/tests/unit/fitting/minimizers/test_minimizer_bumps.py b/tests/unit/fitting/minimizers/test_minimizer_bumps.py index 28d4a049..beb6301c 100644 --- a/tests/unit/fitting/minimizers/test_minimizer_bumps.py +++ b/tests/unit/fitting/minimizers/test_minimizer_bumps.py @@ -41,6 +41,13 @@ def test_all_methods(self, minimizer: Bumps) -> None: # When Then Expect assert minimizer.all_methods() == ['amoeba', 'de', 'dream', 'newton', 'lm'] + def test_all_methods_returns_a_copy(self, minimizer: Bumps) -> None: + """Callers must not be able to mutate the module-level list in place.""" + methods = minimizer.all_methods() + methods.append('tampered') + + assert 'tampered' not in minimizer.all_methods() + def test_supported_methods(self, minimizer: Bumps) -> None: # When Then Expect assert set(minimizer.supported_methods()) == set(['newton', 'lm', 'amoeba']) @@ -65,17 +72,17 @@ def test_fit(self, minimizer: Bumps, monkeypatch) -> None: # Prepare a mock parameter with .name = 'pmock_parm_1' mock_bumps_param = MagicMock() mock_bumps_param.name = 'pmock_parm_1' - # Patch FitProblem to have _parameters attribute as expected - mock_FitProblem_instance = MagicMock() - mock_FitProblem_instance._parameters = [mock_bumps_param] - mock_FitProblem = MagicMock(return_value=mock_FitProblem_instance) + # A mock problem with _parameters, plus the Curve model returned + # directly by the helper (never via the deprecated problem.fitness) + mock_model = MagicMock() + mock_problem = MagicMock() + mock_problem._parameters = [mock_bumps_param] + mock_counter = MagicMock() + mock_build = MagicMock(return_value=(mock_problem, mock_counter, mock_model)) monkeypatch.setattr( - easyscience.fitting.minimizers.minimizer_bumps, 'FitProblem', mock_FitProblem + easyscience.fitting.minimizers.minimizer_bumps, 'build_curve_problem', mock_build ) - mock_model = MagicMock() - mock_model_function = MagicMock(return_value=mock_model) - minimizer._make_model = MagicMock(return_value=mock_model_function) minimizer._gen_fit_results = MagicMock(return_value='gen_fit_results') cached_par = MagicMock() @@ -84,11 +91,10 @@ def test_fit(self, minimizer: Bumps, monkeypatch) -> None: minimizer._cached_pars = cached_pars minimizer._cached_pars_vals = {'mock_parm_1': (1, 0.0)} - # Patch _set_parameter_fit_result - def fake_set_parameter_fit_result(fit_result, stack_status, par_list): - for index, name in enumerate([par.name for par in par_list]): - dict_name = name[len('p') :] - minimizer._cached_pars[dict_name].value = fit_result.x[index] + # Patch _set_parameter_fit_result. It now receives prefix-stripped names. + def fake_set_parameter_fit_result(fit_result, stack_status, par_names): + for index, name in enumerate(par_names): + minimizer._cached_pars[name].value = fit_result.x[index] minimizer._set_parameter_fit_result = fake_set_parameter_fit_result @@ -104,7 +110,17 @@ def fake_set_parameter_fit_result(fit_result, stack_status, par_list): mock_FitDriver.assert_called_once() mock_driver_instance.clip.assert_called_once() mock_driver_instance.fit.assert_called_once() - minimizer._make_model.assert_called_once_with(parameters=None) + # The problem is built via the shared helper and its Curve is cached + mock_build.assert_called_once() + build_args = mock_build.call_args + assert build_args.args[0] is minimizer + assert np.array_equal(build_args.args[1], np.asarray(1.0)) + assert np.array_equal(build_args.args[2], np.asarray(2.0)) + assert np.array_equal(build_args.args[3], np.asarray(1)) + assert build_args.kwargs == {'parameters': None} + assert minimizer._eval_counter is mock_counter + assert minimizer._cached_model is mock_model + assert mock_FitDriver.call_args.kwargs['problem'] is mock_problem # _gen_fit_results is called with the OptimizeResult built from driver.fit() minimizer._gen_fit_results.assert_called_once() passed_result = minimizer._gen_fit_results.call_args.args[0] @@ -114,8 +130,6 @@ def fake_set_parameter_fit_result(fit_result, stack_status, par_list): 'max_evaluations': None, 'tolerance': None, } - mock_model_function.assert_called_once_with(1.0, 2.0, 1) - mock_FitProblem.assert_called_once_with(mock_model) @pytest.mark.parametrize( 'weights', @@ -134,34 +148,6 @@ def test_fit_weight_exceptions(self, minimizer: Bumps, weights) -> None: with pytest.raises(ValueError): minimizer.fit(x=np.array([1, 2, 3]), y=np.array([1, 2, 3]), weights=weights) - def test_make_model(self, minimizer: Bumps, monkeypatch) -> None: - # When - mock_fit_function = MagicMock(return_value=np.array([11, 22])) - minimizer._generate_fit_function = MagicMock(return_value=mock_fit_function) - - mock_parm_1 = MagicMock() - mock_parm_1.unique_name = 'mock_parm_1' - minimizer.convert_to_par_object = MagicMock(return_value='converted_parm_1') - - mock_Curve = MagicMock(return_value='curve') - monkeypatch.setattr(easyscience.fitting.minimizers.minimizer_bumps, 'Curve', mock_Curve) - - # Then - model = minimizer._make_model(parameters=[mock_parm_1]) - curve_for_model = model( - x=np.array([1, 2]), y=np.array([10, 20]), weights=np.array([100, 200]) - ) - wrapped_fit_function = mock_Curve.call_args[0][0] - wrapped_fit_function(np.array([1, 2]), pmock_parm_1=3) - - # Expect - minimizer._generate_fit_function.assert_called_once_with() - assert minimizer._eval_counter is wrapped_fit_function - assert minimizer._eval_counter.count == 1 - assert all(mock_Curve.call_args[0][1] == np.array([1, 2])) - assert all(mock_Curve.call_args[0][2] == np.array([10, 20])) - assert curve_for_model == 'curve' - def test_set_parameter_fit_result_no_stack_status(self, minimizer: Bumps): # When minimizer._cached_pars = { @@ -179,15 +165,8 @@ def test_set_parameter_fit_result_no_stack_status(self, minimizer: Bumps): mock_fit_result.x = np.array([1.0, 2.0]) mock_fit_result.dx = np.array([0.1, 0.2]) - # The new argument: par_list (list of mock parameters) - mock_par_a = MagicMock() - mock_par_a.name = 'pa' - mock_par_b = MagicMock() - mock_par_b.name = 'pb' - par_list = [mock_par_a, mock_par_b] - - # Then - minimizer._set_parameter_fit_result(mock_fit_result, False, par_list) + # Then - names arrive already stripped of the minimizer prefix + minimizer._set_parameter_fit_result(mock_fit_result, False, ['a', 'b']) # Expect assert minimizer._cached_pars['a'].value == 1.0 @@ -195,6 +174,20 @@ def test_set_parameter_fit_result_no_stack_status(self, minimizer: Bumps): assert minimizer._cached_pars['b'].value == 2.0 assert minimizer._cached_pars['b'].error == 0.2 + def test_set_parameter_fit_result_without_stderr(self, minimizer: Bumps): + """Fitters that cannot produce a covariance hand back ``dx=None``; + those parameters get a zero error instead of raising.""" + minimizer._cached_pars = {'a': MagicMock()} + + mock_fit_result = MagicMock() + mock_fit_result.x = np.array([1.0]) + mock_fit_result.dx = None + + minimizer._set_parameter_fit_result(mock_fit_result, False, ['a']) + + assert minimizer._cached_pars['a'].value == 1.0 + assert minimizer._cached_pars['a'].error == 0.0 + def test_gen_fit_results( self, minimizer: Bumps, monkeypatch, caplog: 'pytest.LogCaptureFixture' ): @@ -302,6 +295,51 @@ def test_gen_fit_results_max_evaluations_boundary( assert mock_domain_fit_results.success is expected_success + def test_gen_fit_results_applies_extra_kwargs(self, minimizer: Bumps) -> None: + """Extra kwargs land on a real FitResults. Guarding the copy on the + current value instead of `hasattr` would drop every one of them, since + all FitResults fields start out falsy.""" + mock_cached_model = MagicMock() + mock_cached_model.x = np.array([1.0]) + mock_cached_model.y = np.array([2.0]) + mock_cached_model.dy = np.array([1.0]) + mock_cached_model.pars = {'ppar_1': 0} + minimizer._cached_model = mock_cached_model + minimizer._cached_pars = {'par_1': MagicMock(value=1.0)} + minimizer._p_0 = {} + minimizer._eval_counter = None + minimizer.evaluate = MagicMock(return_value=np.array([2.0])) + + mock_fit_result = MagicMock() + mock_fit_result.success = True + mock_fit_result.nit = 1 + + results = minimizer._gen_fit_results(mock_fit_result, x_matrices='copied') + + assert results.x_matrices == 'copied' + + def test_gen_fit_results_propagates_failure_message(self, minimizer: Bumps) -> None: + mock_cached_model = MagicMock() + mock_cached_model.x = np.array([1.0]) + mock_cached_model.y = np.array([2.0]) + mock_cached_model.dy = np.array([1.0]) + mock_cached_model.pars = {'ppar_1': 0} + minimizer._cached_model = mock_cached_model + minimizer._cached_pars = {'par_1': MagicMock(value=1.0)} + minimizer._p_0 = {} + minimizer._eval_counter = None + minimizer.evaluate = MagicMock(return_value=np.array([2.0])) + + mock_fit_result = MagicMock() + mock_fit_result.success = False + mock_fit_result.nit = 1 + mock_fit_result.message = 'Fit aborted before convergence' + + results = minimizer._gen_fit_results(mock_fit_result) + + assert results.success is False + assert results.message == 'Fit aborted before convergence' + def test_resolve_fitclass_valid(self, minimizer: Bumps) -> None: # When Then fitclass = Bumps._resolve_fitclass('lm') @@ -334,16 +372,14 @@ def test_fit_progress_callback(self, minimizer: Bumps, monkeypatch) -> None: mock_bumps_param = MagicMock() mock_bumps_param.name = 'pmock_parm_1' - mock_FitProblem_instance = MagicMock() - mock_FitProblem_instance._parameters = [mock_bumps_param] - mock_FitProblem = MagicMock(return_value=mock_FitProblem_instance) + mock_problem = MagicMock() + mock_problem._parameters = [mock_bumps_param] monkeypatch.setattr( - easyscience.fitting.minimizers.minimizer_bumps, 'FitProblem', mock_FitProblem + easyscience.fitting.minimizers.minimizer_bumps, + 'build_curve_problem', + MagicMock(return_value=(mock_problem, MagicMock(), MagicMock())), ) - mock_model = MagicMock() - mock_model_function = MagicMock(return_value=mock_model) - minimizer._make_model = MagicMock(return_value=mock_model_function) minimizer._set_parameter_fit_result = MagicMock() minimizer._gen_fit_results = MagicMock(return_value='gen_fit_results') @@ -363,7 +399,7 @@ def test_fit_progress_callback(self, minimizer: Bumps, monkeypatch) -> None: monitors = driver_call_kwargs.kwargs.get('monitors', driver_call_kwargs[1].get('monitors')) assert len(monitors) == 1 assert isinstance(monitors[0], BumpsProgressMonitor) - assert monitors[0]._problem is mock_FitProblem_instance + assert monitors[0]._problem is mock_problem assert monitors[0]._callback is progress_callback assert monitors[0]._payload_builder == minimizer._build_progress_payload @@ -394,12 +430,21 @@ def test_fit_uses_supplied_model_and_optional_kwargs( MagicMock(return_value=mock_problem), ) - minimizer._make_model = MagicMock() + mock_build = MagicMock() + monkeypatch.setattr( + easyscience.fitting.minimizers.minimizer_bumps, 'build_curve_problem', mock_build + ) minimizer._gen_fit_results = MagicMock(return_value='gen_fit_results') minimizer._resolve_fitclass = MagicMock(return_value=MagicMock(id='amoeba')) minimizer._set_parameter_fit_result = MagicMock() - minimizer._cached_pars = {'mock_parm_1': MagicMock(value=1.0)} - minimizer._cached_pars_vals = {'mock_parm_1': (1.0, 0.0)} + + # A supplied model bypasses build_curve_problem, so fit() must populate the + # parameter cache itself from the bound object rather than leaving it empty. + object_parameter = MagicMock(unique_name='mock_parm_1') + object_parameter.value = 1.0 + object_parameter.error = 0.0 + minimizer._object = MagicMock() + minimizer._object.get_fit_parameters = MagicMock(return_value=[object_parameter]) supplied_model = MagicMock() minimizer_kwargs = {'existing_option': 'minimizer'} @@ -417,7 +462,7 @@ def test_fit_uses_supplied_model_and_optional_kwargs( ) assert result == 'gen_fit_results' - minimizer._make_model.assert_not_called() + mock_build.assert_not_called() fit_driver_kwargs = mock_FitDriver.call_args.kwargs assert fit_driver_kwargs['problem'] is mock_problem assert fit_driver_kwargs['existing_option'] == 'minimizer' @@ -426,6 +471,53 @@ def test_fit_uses_supplied_model_and_optional_kwargs( assert fit_driver_kwargs['xtol'] == 0.25 assert fit_driver_kwargs['steps'] == 7 mock_driver_instance.fit.assert_called_once() + # The cache and the starting-point snapshot are built from the bound object + assert minimizer._cached_pars == {'mock_parm_1': object_parameter} + assert minimizer._p_0 == {'pmock_parm_1': 1.0} + + def test_fit_with_supplied_model_resets_eval_counter( + self, minimizer: Bumps, monkeypatch + ) -> None: + """A supplied model installs no EvalCounter, so a counter left over + from a previous fit must not be reported as this fit's count.""" + from easyscience import global_object + + global_object.stack.enabled = False + + mock_driver_instance = MagicMock() + mock_driver_instance.fit = MagicMock(return_value=(np.array([3.0]), 0.0)) + mock_driver_instance.stderr = MagicMock(return_value=np.array([0.1])) + mock_driver_instance.monitor_runner.history.step = [0] + monkeypatch.setattr( + easyscience.fitting.minimizers.minimizer_bumps, + 'FitDriver', + MagicMock(return_value=mock_driver_instance), + ) + mock_problem = MagicMock() + mock_problem._parameters = [] + monkeypatch.setattr( + easyscience.fitting.minimizers.minimizer_bumps, + 'FitProblem', + MagicMock(return_value=mock_problem), + ) + + minimizer._gen_fit_results = MagicMock(return_value='gen_fit_results') + minimizer._resolve_fitclass = MagicMock(return_value=MagicMock(id='amoeba')) + minimizer._set_parameter_fit_result = MagicMock() + minimizer._object = MagicMock() + minimizer._object.get_fit_parameters = MagicMock(return_value=[]) + + # Stale counter from an earlier fit + minimizer._eval_counter = MagicMock(count=999) + + minimizer.fit( + x=np.array([1.0]), + y=np.array([2.0]), + weights=np.array([1.0]), + model=MagicMock(), + ) + + assert minimizer._eval_counter is None def test_fit_rejects_non_callable_progress_callback( self, minimizer: Bumps, monkeypatch @@ -522,19 +614,6 @@ def test_build_progress_payload_reduced_chi2_positive_dof(self, minimizer: Bumps ((), {'nllf': 5.0, 'norm': True}), ] - def test_current_parameter_snapshot(self, minimizer: Bumps) -> None: - # When - mock_problem = MagicMock() - mock_problem.labels.return_value = ['palpha', 'pbeta'] - - point = np.array([1.5, 2.5]) - - # Then - snapshot = minimizer._current_parameter_snapshot(mock_problem, point) - - # Expect - assert snapshot == {'alpha': 1.5, 'beta': 2.5} - @pytest.mark.parametrize('par_list', [None, [MagicMock(unique_name='alpha')]]) def test_convert_to_pars_obj_optional_parameter_list( self, minimizer: Bumps, par_list, monkeypatch @@ -557,25 +636,6 @@ def test_convert_to_pars_obj_optional_parameter_list( else: minimizer._object.get_fit_parameters.assert_not_called() - def test_make_model_without_parameters_uses_cached_parameters( - self, minimizer: Bumps, monkeypatch - ) -> None: - minimizer._generate_fit_function = MagicMock( - return_value=MagicMock(return_value=np.array([2.0])) - ) - minimizer._cached_pars = {'alpha': MagicMock(value=1.0)} - minimizer.convert_to_par_object = MagicMock(return_value='converted-alpha') - - mock_curve = MagicMock(return_value='curve') - monkeypatch.setattr(easyscience.fitting.minimizers.minimizer_bumps, 'Curve', mock_curve) - - model = minimizer._make_model() - curve = model(np.array([1.0]), np.array([2.0]), np.array([3.0])) - - assert curve == 'curve' - minimizer.convert_to_par_object.assert_called_once_with(minimizer._cached_pars['alpha']) - assert mock_curve.call_args.kwargs['palpha'] == 'converted-alpha' - def test_bumps_progress_monitor_calls_callback(self, minimizer: Bumps) -> None: # When callback = MagicMock(return_value=True) @@ -622,16 +682,13 @@ def test_fit_exception_restores_values(self, minimizer: Bumps, monkeypatch) -> N easyscience.fitting.minimizers.minimizer_bumps, 'FitDriver', mock_FitDriver ) - mock_FitProblem_instance = MagicMock() - mock_FitProblem_instance._parameters = [] - mock_FitProblem = MagicMock(return_value=mock_FitProblem_instance) + mock_problem = MagicMock() + mock_problem._parameters = [] monkeypatch.setattr( - easyscience.fitting.minimizers.minimizer_bumps, 'FitProblem', mock_FitProblem + easyscience.fitting.minimizers.minimizer_bumps, + 'build_curve_problem', + MagicMock(return_value=(mock_problem, MagicMock(), MagicMock())), ) - - mock_model = MagicMock() - mock_model_function = MagicMock(return_value=mock_model) - minimizer._make_model = MagicMock(return_value=mock_model_function) minimizer._resolve_fitclass = MagicMock(return_value=MagicMock(id='amoeba')) # Then Expect @@ -683,532 +740,214 @@ def test_gen_fit_results_uses_nit_for_budget_check( # =================================================================== -# Bumps.mcmc_sample() — Bayesian DREAM sampling +# fit() — tolerance / budget defaults are reported, never forced # =================================================================== -class TestBumpsSample: - """Tests for the ``Bumps.mcmc_sample()`` method and its helpers.""" - - # Sentinel value to signal "set fitter.state = None" in _setup_driver_mock - ABORT = object() +class TestFitToleranceAndBudgetDefaults: + """BUMPS pairs an independent ftol/xtol default per fitter. Resolving + them for reporting must not push a single collapsed value back into the + fitter, which would silently tighten its convergence criteria.""" @pytest.fixture def minimizer(self) -> Bumps: return Bumps( obj='obj', fit_function='fit_function', - minimizer_enum=MagicMock(package='bumps', method='amoeba'), - ) - - @pytest.fixture(autouse=True) - def _mock_bumps_internals(self, monkeypatch): - """Prevent sample() from constructing real BUMPS objects. - - ``sample()`` imports ``DreamFit`` and ``FitProblem`` from the real - ``bumps`` package internally, which would try to build real model - objects. We redirect those to mocks and also mock ``FitDriver`` - (which *is* a module-level import) so the whole flow stays under - test control. - - Also mock ``_make_model`` on the class so that the ``minimizer`` - fixture (which uses ``obj='obj'``) doesn't fail inside ``sample()``. - """ - import bumps.fitters - import bumps.names - - monkeypatch.setattr(bumps.fitters, 'DreamFit', MagicMock()) - monkeypatch.setattr(bumps.names, 'FitProblem', MagicMock(return_value=MagicMock())) - monkeypatch.setattr( - Bumps, '_make_model', MagicMock(return_value=MagicMock(return_value=MagicMock())) + minimizer_enum=MagicMock(package='bumps', method='newton'), ) - def _setup_driver_mock( - self, monkeypatch, fitter_state_value=None, fit_result=None, fit_side_effect=None - ): - """Helper to create a mocked FitDriver with configurable behavior. - - :param fitter_state_value: If ``None``, ``driver.fitter.state`` will be - a regular MagicMock (non-None). Pass ``ABORT`` to set it to ``None`` - and simulate user abort. - """ + @staticmethod + def _patch_driver_and_problem(minimizer: Bumps, monkeypatch) -> MagicMock: from easyscience import global_object global_object.stack.enabled = False mock_driver = MagicMock() - mock_driver.clip = MagicMock() - - if fit_side_effect is not None: - mock_driver.fit.side_effect = fit_side_effect - else: - mock_driver.fit.return_value = fit_result or (np.array([1.0]), 0.0) - + mock_driver.fit = MagicMock(return_value=(np.array([42.0]), 0.0)) mock_driver.stderr = MagicMock(return_value=np.array([0.1])) - - if fitter_state_value is TestBumpsSample.ABORT: - mock_driver.fitter.state = None - else: - mock_state = MagicMock() - mock_state.Nvar = 1 - mock_state.Npop = 5 - mock_state.labels = ['p_param_0'] - mock_draw = MagicMock() - mock_draw.points = np.array([[1.0]]) - mock_draw.logp = np.array([0.5]) - mock_state.draw.return_value = mock_draw - mock_driver.fitter.state = mock_state - + mock_driver.monitor_runner.history.step = [0] mock_FitDriver = MagicMock(return_value=mock_driver) monkeypatch.setattr( easyscience.fitting.minimizers.minimizer_bumps, 'FitDriver', mock_FitDriver ) - return mock_FitDriver, mock_driver - @pytest.mark.parametrize( - 'kwargs, match', - [ - ({'samples': 0}, 'samples must be a positive integer'), - ({'samples': -1}, 'samples must be a positive integer'), - ({'burn': -1}, 'burn must be a non-negative integer'), - ({'thin': 0}, 'thin must be a positive integer'), - ], - ) - def test_sample_invalid_args(self, minimizer: Bumps, kwargs, match) -> None: - """Invalid samples/burn/thin values raise ValueError before any sampling. - - This is the single source of truth for these checks — the higher-level - ``Sampler`` relies on it. - """ - with pytest.raises(ValueError, match=match): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=kwargs.get('samples', 10), - burn=kwargs.get('burn', 0), - thin=kwargs.get('thin', 1), - ) - - @pytest.mark.parametrize( - 'overrides, match', - [ - ({'y': np.array([0.1])}, 'x and y must have the same shape'), - ({'x': np.array([1.0, np.nan])}, 'x cannot contain NaN'), - ({'y': np.array([0.1, np.inf])}, 'y cannot contain NaN'), - ({'weights': np.array([1.0])}, 'Weights must have the same shape'), - ({'weights': np.array([1.0, np.nan])}, 'Weights cannot be NaN'), - ({'weights': np.array([1.0, 0.0])}, 'Weights must be strictly positive'), - ], - ) - def test_sample_invalid_data(self, minimizer: Bumps, overrides, match) -> None: - """Shape mismatches and non-finite/non-positive data raise ValueError - before any sampling.""" - data = { - 'x': np.array([1.0, 2.0]), - 'y': np.array([0.1, 0.2]), - 'weights': np.array([1.0, 1.0]), - } - data.update(overrides) - with pytest.raises(ValueError, match=match): - minimizer.mcmc_sample(**data, samples=10, burn=0, thin=1) - - def test_sample_basic(self, minimizer: Bumps, monkeypatch) -> None: - """Verify that mcmc_sample() returns a dict with expected keys.""" - mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - - result = minimizer.mcmc_sample( - x=np.array([1.0, 2.0]), - y=np.array([0.1, 0.2]), - weights=np.array([1.0, 1.0]), - samples=100, - burn=20, - thin=2, - population=5, - ) - - assert isinstance(result, dict) - assert 'draws' in result - assert 'param_names' in result - assert 'internal_bumps_object' in result - assert 'logp' in result - mock_FitDriver.assert_called_once() - - def test_sample_with_progress_callback(self, minimizer: Bumps, monkeypatch) -> None: - """Verify progress callback is wired up as a monitor.""" - mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - progress_callback = MagicMock() - - result = minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=5, - thin=1, - progress_callback=progress_callback, + mock_problem = MagicMock() + mock_problem._parameters = [] + monkeypatch.setattr( + easyscience.fitting.minimizers.minimizer_bumps, + 'build_curve_problem', + MagicMock(return_value=(mock_problem, MagicMock(count=3), MagicMock())), ) - assert result is not None - call_kwargs = mock_FitDriver.call_args.kwargs - assert 'monitors' in call_kwargs - assert len(call_kwargs['monitors']) == 1 - assert isinstance(call_kwargs['monitors'][0], BumpsProgressMonitor) - - def test_sample_aborted_by_user_raises_fit_error(self, minimizer: Bumps, monkeypatch) -> None: - """Verify that sampling abortion raises FitError.""" - self._setup_driver_mock(monkeypatch, fitter_state_value=TestBumpsSample.ABORT) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - - with pytest.raises(FitError, match='Sampling aborted by user'): - minimizer.mcmc_sample(x=np.array([1.0]), y=np.array([0.1]), weights=np.array([1.0])) + minimizer._gen_fit_results = MagicMock(return_value='result') + minimizer._set_parameter_fit_result = MagicMock() + minimizer._cached_pars = {} + minimizer._cached_pars_vals = {} + return mock_FitDriver - def test_sample_driver_exception_restores_parameters( + def test_tolerance_none_does_not_override_fitter_defaults( self, minimizer: Bumps, monkeypatch ) -> None: - """Verify that a driver exception during sampling restores parameter values.""" - self._setup_driver_mock(monkeypatch, fit_side_effect=RuntimeError('driver failed')) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - minimizer._restore_parameter_values = MagicMock() + mock_FitDriver = self._patch_driver_and_problem(minimizer, monkeypatch) - with pytest.raises(RuntimeError, match='driver failed'): - minimizer.mcmc_sample(x=np.array([1.0]), y=np.array([0.1]), weights=np.array([1.0])) + minimizer.fit(x=np.array([1.0]), y=np.array([2.0]), weights=np.array([1.0])) - minimizer._restore_parameter_values.assert_called_once() + # The real 'newton' settings are ftol=1e-6 / xtol=1e-12. Neither may be + # forwarded, or BUMPS would run against a tolerance the caller never asked for. + driver_kwargs = mock_FitDriver.call_args.kwargs + assert 'ftol' not in driver_kwargs + assert 'xtol' not in driver_kwargs + assert 'steps' not in driver_kwargs - def test_sample_population_param(self, minimizer: Bumps, monkeypatch) -> None: - """population kwarg is forwarded to DREAM as pop.""" - mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) + # ...but the resolved defaults are still reported for the budget check. + gen_kwargs = minimizer._gen_fit_results.call_args.kwargs + assert gen_kwargs['tolerance'] == 1e-12 # min(ftol, xtol) + assert gen_kwargs['max_evaluations'] == 3000 # 'newton' default steps - minimizer.mcmc_sample( + def test_explicit_tolerance_is_forwarded(self, minimizer: Bumps, monkeypatch) -> None: + mock_FitDriver = self._patch_driver_and_problem(minimizer, monkeypatch) + + minimizer.fit( x=np.array([1.0]), - y=np.array([0.1]), + y=np.array([2.0]), weights=np.array([1.0]), - samples=10, - burn=0, - thin=1, - population=7, + tolerance=1e-3, + max_evaluations=11, ) - call_kwargs = mock_FitDriver.call_args.kwargs - assert call_kwargs['pop'] == 7 + driver_kwargs = mock_FitDriver.call_args.kwargs + assert driver_kwargs['ftol'] == 1e-3 + assert driver_kwargs['xtol'] == 1e-3 + assert driver_kwargs['steps'] == 11 - def test_sample_sampler_kwargs_forwarded(self, minimizer: Bumps, monkeypatch) -> None: - """sampler_kwargs entries are merged into the DREAM kwargs.""" - mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) + def test_minimizer_kwargs_is_not_mutated(self, minimizer: Bumps, monkeypatch) -> None: + self._patch_driver_and_problem(minimizer, monkeypatch) - minimizer.mcmc_sample( + minimizer_kwargs = {'existing': 'value'} + minimizer.fit( x=np.array([1.0]), - y=np.array([0.1]), + y=np.array([2.0]), weights=np.array([1.0]), - samples=10, - burn=0, - thin=1, - sampler_kwargs={'trim': False}, + tolerance=1e-3, + max_evaluations=11, + minimizer_kwargs=minimizer_kwargs, + engine_kwargs={'engine': 'option'}, ) - assert mock_FitDriver.call_args.kwargs['trim'] is False - - def test_sample_rejects_non_callable_callback(self, minimizer: Bumps, monkeypatch) -> None: - import bumps.names + # The caller's mapping is untouched, so reusing it cannot leak settings + # from one fit into the next. + assert minimizer_kwargs == {'existing': 'value'} - monkeypatch.setattr(bumps.names, 'FitProblem', MagicMock(return_value=MagicMock())) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - with pytest.raises(ValueError, match='progress_callback must be callable'): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=5, - thin=1, - progress_callback='not-callable', - ) - - # --- Resume-state tests ------------------------------------------------- - - def _make_resume_state_mock(self, *, nvar=2, npop=10, labels=None): - """Build a mock MCMCDraw for resume tests. - - BUMPS labels follow the pattern ``'p'`` (the - ``MINIMIZER_PARAMETER_PREFIX`` concatenated with the unique name), - e.g. ``'pFilm_thickness'``. - - :param nvar: Number of parameters. - :param npop: Population size. - :param labels: Parameter labels (defaults to ``['pa', 'pb']`` - which strip to ``['a', 'b']``). - """ - if labels is None: - labels = ['pa', 'pb'] - mock_state = MagicMock() - mock_state.Nvar = nvar - mock_state.Npop = npop - mock_state.labels = labels - mock_draw = MagicMock() - mock_draw.points = np.ones((20, nvar)) - mock_draw.logp = np.ones(20) - mock_state.draw.return_value = mock_draw - return mock_state - - def _make_problem_with_parameters(self, param_names): - """Build a mock FitProblem whose ``_parameters`` yields the given names.""" - params = [] - for name in param_names: - p = MagicMock() - p.name = 'p' + name - params.append(p) - mock_problem = MagicMock() - mock_problem._parameters = params - return mock_problem - - def test_sample_resume_state(self, minimizer: Bumps, monkeypatch) -> None: - """Verify resume_state is forwarded to driver.fit().""" - mock_FitDriver, mock_driver = self._setup_driver_mock(monkeypatch) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - resume_state = self._make_resume_state_mock() - - import bumps.names +# =================================================================== +# fit() — unsuccessful and aborted outcomes +# =================================================================== - monkeypatch.setattr( - bumps.names, - 'FitProblem', - MagicMock(return_value=self._make_problem_with_parameters(['a', 'b'])), - ) - result = minimizer.mcmc_sample( - x=np.array([1.0, 2.0]), - y=np.array([0.1, 0.2]), - weights=np.array([1.0, 1.0]), - samples=10, - burn=0, - thin=1, - resume_state=resume_state, +class TestFitUnsuccessfulOutcomes: + @pytest.fixture + def minimizer(self) -> Bumps: + return Bumps( + obj='obj', + fit_function='fit_function', + minimizer_enum=MagicMock(package='bumps', method='amoeba'), ) - assert result is not None - # Verify a fit_state (defensive copy of resume_state) was passed to driver.fit() - call_kwargs = mock_driver.fit.call_args.kwargs - assert call_kwargs.get('fit_state') is not None - assert call_kwargs['fit_state'] is not resume_state + @staticmethod + def _patch(minimizer: Bumps, monkeypatch, driver_result, history_step=None) -> MagicMock: + from easyscience import global_object - def test_sample_resume_param_mismatch_raises(self, minimizer: Bumps, monkeypatch) -> None: - """Parameter count mismatch raises ValueError before driver.fit().""" - import bumps.names + global_object.stack.enabled = False + mock_driver = MagicMock() + mock_driver.fit = MagicMock(return_value=driver_result) + mock_driver.stderr = MagicMock(return_value=np.array([0.1])) + mock_driver.monitor_runner.history.step = [] if history_step is None else history_step monkeypatch.setattr( - bumps.names, - 'FitProblem', - MagicMock(return_value=self._make_problem_with_parameters(['a'])), + easyscience.fitting.minimizers.minimizer_bumps, + 'FitDriver', + MagicMock(return_value=mock_driver), ) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - # resume_state has 2 params, model has 1 - resume_state = self._make_resume_state_mock(nvar=2) - - with pytest.raises(ValueError, match='resume_state has 2 parameters'): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=0, - thin=1, - resume_state=resume_state, - ) - - def test_sample_resume_param_name_mismatch_raises(self, minimizer: Bumps, monkeypatch) -> None: - """Parameter name/order mismatch raises ValueError before driver.fit().""" - import bumps.names + mock_problem = MagicMock() + mock_problem._parameters = [] monkeypatch.setattr( - bumps.names, - 'FitProblem', - MagicMock(return_value=self._make_problem_with_parameters(['a', 'b'])), + easyscience.fitting.minimizers.minimizer_bumps, + 'build_curve_problem', + MagicMock(return_value=(mock_problem, MagicMock(count=3), MagicMock())), ) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - # resume_state has labels ['px', 'py'] → stripped to ['x', 'y'] - # Current model has params ['pa', 'pb'] → stripped to ['a', 'b'] - # → mismatch - resume_state = self._make_resume_state_mock(nvar=2, labels=['px', 'py']) - - with pytest.raises(ValueError, match='Parameter names/order mismatch'): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=0, - thin=1, - resume_state=resume_state, - ) - def test_sample_resume_population_mismatch_raises(self, minimizer: Bumps, monkeypatch) -> None: - """Explicit population differing from state.Npop raises ValueError.""" - import bumps.names - - monkeypatch.setattr( - bumps.names, - 'FitProblem', - MagicMock(return_value=self._make_problem_with_parameters(['a', 'b'])), - ) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - resume_state = self._make_resume_state_mock(nvar=2, npop=10) + minimizer._gen_fit_results = MagicMock(return_value='result') + minimizer._set_parameter_fit_result = MagicMock() + minimizer._resolve_fitclass = MagicMock(return_value=MagicMock(id='amoeba')) + minimizer._cached_pars = {} + minimizer._cached_pars_vals = {} + return mock_driver + + def test_no_solution_is_reported_not_raised(self, minimizer: Bumps, monkeypatch) -> None: + """BUMPS returns x=None for a failed optimization (e.g. LM landing on + non-finite values). That is a non-converged fit, not an exception.""" + self._patch(minimizer, monkeypatch, driver_result=(None, None), history_step=[4]) + minimizer._restore_parameter_values = MagicMock() - with pytest.raises(ValueError, match='would produce'): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=0, - thin=1, - population=3, # ceil(3*2)=6 ≠ 10 - resume_state=resume_state, - ) + result = minimizer.fit(x=np.array([1.0]), y=np.array([2.0]), weights=np.array([1.0])) - def test_sample_resume_forces_burn_to_zero( - self, minimizer: Bumps, monkeypatch, caplog: 'pytest.LogCaptureFixture' - ) -> None: - """burn>0 with resume_state warns and is forced to 0.""" - mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) - import bumps.names + assert result == 'result' + passed = minimizer._gen_fit_results.call_args.args[0] + assert passed.success is False + assert passed.x is None + assert passed.dx is None # stderr() needs a solution to expand around + assert 'did not converge' in passed.message + # Parameters are rolled back and never written from a missing solution + minimizer._restore_parameter_values.assert_called_once() + minimizer._set_parameter_fit_result.assert_not_called() - monkeypatch.setattr( - bumps.names, - 'FitProblem', - MagicMock(return_value=self._make_problem_with_parameters(['a', 'b'])), + def test_abort_is_reported_as_unsuccessful(self, minimizer: Bumps, monkeypatch) -> None: + self._patch( + minimizer, monkeypatch, driver_result=(np.array([42.0]), 0.0), history_step=[2] ) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - resume_state = self._make_resume_state_mock() - - with caplog.at_level(logging.WARNING, logger='easyscience.fitting.bumps'): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=5, - thin=1, - resume_state=resume_state, - ) - assert 'ignored on resume' in caplog.text - # burn must be forced to 0 in the kwargs passed to BUMPS - assert mock_FitDriver.call_args.kwargs['burn'] == 0 - - def test_sample_resume_unlabeled_state_warns_and_uses_absolute_pop( - self, minimizer: Bumps, monkeypatch, caplog: 'pytest.LogCaptureFixture' - ) -> None: - """A state reloaded from disk carries default labels ('P0', ...), so - name validation is skipped with a warning, and the saved population is - reproduced as a negative pop (BUMPS' absolute-chain-count convention).""" - mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) - import bumps.names - - monkeypatch.setattr( - bumps.names, - 'FitProblem', - MagicMock(return_value=self._make_problem_with_parameters(['a', 'b'])), + result = minimizer.fit( + x=np.array([1.0]), + y=np.array([2.0]), + weights=np.array([1.0]), + abort_test=lambda: True, ) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) - resume_state = self._make_resume_state_mock(nvar=2, npop=10, labels=['P0', 'P1']) - with caplog.at_level(logging.WARNING, logger='easyscience.fitting.bumps'): - minimizer.mcmc_sample( - x=np.array([1.0]), - y=np.array([0.1]), - weights=np.array([1.0]), - samples=10, - burn=0, - thin=1, - resume_state=resume_state, - ) + assert result == 'result' + passed = minimizer._gen_fit_results.call_args.args[0] + assert passed.success is False + assert passed.message == 'Fit aborted before convergence' + # The best point reached before the abort is still applied + minimizer._set_parameter_fit_result.assert_called_once() - assert 'does not carry parameter names' in caplog.text - assert mock_FitDriver.call_args.kwargs['pop'] == -10 + def test_empty_step_history_does_not_raise(self, minimizer: Bumps, monkeypatch) -> None: + """An abort before the fitter reports its first step leaves the BUMPS + history trace empty; indexing it would raise IndexError.""" + self._patch(minimizer, monkeypatch, driver_result=(np.array([42.0]), 0.0), history_step=[]) + result = minimizer.fit(x=np.array([1.0]), y=np.array([2.0]), weights=np.array([1.0])) -# =================================================================== -# _build_sample_progress_payload -# =================================================================== + assert result == 'result' + assert minimizer._gen_fit_results.call_args.args[0].nit is None - -class TestBuildSampleProgressPayload: - @pytest.fixture - def minimizer(self) -> Bumps: - return Bumps( - obj='obj', - fit_function='fit_function', - minimizer_enum=MagicMock(package='bumps', method='amoeba'), + def test_successful_fit_reports_success(self, minimizer: Bumps, monkeypatch) -> None: + self._patch( + minimizer, monkeypatch, driver_result=(np.array([42.0]), 0.0), history_step=[7] ) - def test_payload_structure_and_sampling_flag(self, minimizer: Bumps) -> None: - b = minimizer - - mock_problem = MagicMock() - mock_problem.chisq.side_effect = [25.0, 12.5] - mock_problem.labels.return_value = ['palpha'] - mock_problem.getp.return_value = np.array([1.0]) - b._cached_pars = {'alpha': MagicMock(value=1.0)} - - payload = b._build_sample_progress_payload(mock_problem, 7, np.array([1.0]), 12.5) - - assert payload['iteration'] == 7 - assert payload['chi2'] == 25.0 - assert payload['reduced_chi2'] == 12.5 - assert payload['parameter_values'] == {'alpha': 1.0} - assert payload['sampling'] is True - assert payload['finished'] is False - assert payload['refresh_plots'] is False - - def test_payload_keys(self, minimizer: Bumps) -> None: - b = minimizer - mock_problem = MagicMock() - mock_problem.chisq.side_effect = [10.0, 5.0] - mock_problem.labels.return_value = ['pa'] - mock_problem.getp.return_value = np.array([5.0]) - b._cached_pars = {'a': MagicMock(value=5.0)} - - payload = b._build_sample_progress_payload(mock_problem, 1, np.array([5.0]), nllf=5.0) - - expected_keys = { - 'iteration', - 'chi2', - 'reduced_chi2', - 'parameter_values', - 'refresh_plots', - 'finished', - 'sampling', - } - assert set(payload.keys()) == expected_keys - - def test_delegates_to_build_progress_payload(self, minimizer: Bumps) -> None: - """_build_sample_progress_payload calls _build_progress_payload and adds sampling.""" - mock_problem = MagicMock() - - # Patch _build_progress_payload to track calls - base_payload = { - 'iteration': 3, - 'chi2': 42.0, - 'reduced_chi2': 21.0, - 'parameter_values': {'x': 7.0}, - 'refresh_plots': False, - 'finished': False, - } - with patch.object( - minimizer, '_build_progress_payload', return_value=base_payload - ) as mock_bpp: - result = minimizer._build_sample_progress_payload( - mock_problem, 3, np.array([7.0]), 21.0 - ) + minimizer.fit( + x=np.array([1.0]), + y=np.array([2.0]), + weights=np.array([1.0]), + abort_test=lambda: False, + ) - mock_bpp.assert_called_once_with(mock_problem, 3, np.array([7.0]), 21.0) - assert result == {**base_payload, 'sampling': True} + passed = minimizer._gen_fit_results.call_args.args[0] + assert passed.success is True + assert passed.message == 'successful termination' + assert passed.nit == 7 # =================================================================== @@ -1239,13 +978,7 @@ def test_stack_status_true_calls_begin_end_macro(self, minimizer: Bumps) -> None mock_fit_result.x = np.array([1.0, 2.0]) mock_fit_result.dx = np.array([0.1, 0.2]) - mock_par_a = MagicMock() - mock_par_a.name = 'pa' - mock_par_b = MagicMock() - mock_par_b.name = 'pb' - par_list = [mock_par_a, mock_par_b] - - minimizer._set_parameter_fit_result(mock_fit_result, True, par_list) + minimizer._set_parameter_fit_result(mock_fit_result, True, ['a', 'b']) assert minimizer._cached_pars['a'].value == 1.0 assert minimizer._cached_pars['a'].error == 0.1 @@ -1319,11 +1052,10 @@ def test_abort_test_passed_to_fit_driver(self, minimizer: Bumps, monkeypatch) -> mock_problem._parameters = [] monkeypatch.setattr( easyscience.fitting.minimizers.minimizer_bumps, - 'FitProblem', - MagicMock(return_value=mock_problem), + 'build_curve_problem', + MagicMock(return_value=(mock_problem, MagicMock(), MagicMock())), ) - minimizer._make_model = MagicMock(return_value=MagicMock(return_value=MagicMock())) minimizer._gen_fit_results = MagicMock(return_value='result') minimizer._resolve_fitclass = MagicMock(return_value=MagicMock(id='amoeba')) minimizer._set_parameter_fit_result = MagicMock() diff --git a/tests/unit/fitting/samplers/__init__.py b/tests/unit/fitting/samplers/__init__.py new file mode 100644 index 00000000..4e798e20 --- /dev/null +++ b/tests/unit/fitting/samplers/__init__.py @@ -0,0 +1,2 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause diff --git a/tests/unit/fitting/samplers/test_sampler_dream.py b/tests/unit/fitting/samplers/test_sampler_dream.py new file mode 100644 index 00000000..732f5cab --- /dev/null +++ b/tests/unit/fitting/samplers/test_sampler_dream.py @@ -0,0 +1,483 @@ +# SPDX-FileCopyrightText: 2026 EasyScience contributors +# SPDX-License-Identifier: BSD-3-Clause +"""Unit tests for ``DreamSampler`` — mirrors +``src/easyscience/fitting/samplers/sampler_dream.py``. + +Ported from the former ``TestBumpsSample`` suite in +``tests/unit/fitting/minimizers/test_minimizer_bumps.py`` when the former +``Bumps.mcmc_sample`` moved here as ``DreamSampler.run`` (easyscience/core#280). +That entry point, and ``Fitter.mcmc_sample``, have since been removed — +``Sampler`` and ``DreamSampler.run`` are the supported APIs. +""" + +import logging +from unittest.mock import MagicMock + +import numpy as np +import pytest + +import easyscience.fitting.samplers.sampler_dream +from easyscience.fitting.engine_base import EngineBase +from easyscience.fitting.minimizers.bumps_utils import BumpsProgressMonitor +from easyscience.fitting.minimizers.utils import FitError +from easyscience.fitting.samplers import DreamSampler + + +class TestDreamSamplerRun: + """Tests for ``DreamSampler.run()`` and its helpers.""" + + # Sentinel value to signal "set fitter.state = None" in _setup_driver_mock + ABORT = object() + + @pytest.fixture + def engine(self) -> DreamSampler: + return DreamSampler(obj='obj', fit_function='fit_function') + + @pytest.fixture(autouse=True) + def _mock_bumps_internals(self, monkeypatch): + """Prevent run() from constructing real BUMPS objects. + + ``run()`` imports ``DreamFit`` from the real ``bumps`` package + internally and builds its problem via ``build_curve_problem``, + which would try to build real model objects. We redirect those + to mocks and also mock ``FitDriver`` (a module-level import) so + the whole flow stays under test control. + """ + import bumps.fitters + + monkeypatch.setattr(bumps.fitters, 'DreamFit', MagicMock()) + self._set_problem(monkeypatch, MagicMock()) + + @staticmethod + def _set_problem(monkeypatch, problem): + """Point ``build_curve_problem`` at a canned (problem, counter, curve) triple.""" + monkeypatch.setattr( + easyscience.fitting.samplers.sampler_dream, + 'build_curve_problem', + MagicMock(return_value=(problem, MagicMock(), MagicMock())), + ) + + def _setup_driver_mock( + self, monkeypatch, fitter_state_value=None, fit_result=None, fit_side_effect=None + ): + """Helper to create a mocked FitDriver with configurable behavior. + + :param fitter_state_value: If ``None``, ``driver.fitter.state`` will be + a regular MagicMock (non-None). Pass ``ABORT`` to set it to ``None`` + and simulate user abort. + """ + from easyscience import global_object + + global_object.stack.enabled = False + + mock_driver = MagicMock() + mock_driver.clip = MagicMock() + + if fit_side_effect is not None: + mock_driver.fit.side_effect = fit_side_effect + else: + mock_driver.fit.return_value = fit_result or (np.array([1.0]), 0.0) + + mock_driver.stderr = MagicMock(return_value=np.array([0.1])) + + if fitter_state_value is TestDreamSamplerRun.ABORT: + mock_driver.fitter.state = None + else: + mock_state = MagicMock() + mock_state.Nvar = 1 + mock_state.Npop = 5 + mock_state.labels = ['p_param_0'] + mock_draw = MagicMock() + mock_draw.points = np.array([[1.0]]) + mock_draw.logp = np.array([0.5]) + mock_state.draw.return_value = mock_draw + mock_driver.fitter.state = mock_state + + mock_FitDriver = MagicMock(return_value=mock_driver) + monkeypatch.setattr( + easyscience.fitting.samplers.sampler_dream, 'FitDriver', mock_FitDriver + ) + return mock_FitDriver, mock_driver + + def test_is_an_engine(self, engine: DreamSampler) -> None: + assert isinstance(engine, EngineBase) + assert engine.package == 'bumps' + + @pytest.mark.parametrize( + 'kwargs, match', + [ + ({'samples': 0}, 'samples must be a positive integer'), + ({'samples': -1}, 'samples must be a positive integer'), + ({'burn': -1}, 'burn must be a non-negative integer'), + ({'thin': 0}, 'thin must be a positive integer'), + ], + ) + def test_run_invalid_args(self, engine: DreamSampler, kwargs, match) -> None: + """Invalid samples/burn/thin values raise ValueError before any sampling. + + This is the single source of truth for these checks — the higher-level + ``Sampler`` relies on it. + """ + with pytest.raises(ValueError, match=match): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=kwargs.get('samples', 10), + burn=kwargs.get('burn', 0), + thin=kwargs.get('thin', 1), + ) + + @pytest.mark.parametrize( + 'overrides, match', + [ + ({'y': np.array([0.1])}, 'x and y must have the same shape'), + ({'x': np.array([1.0, np.nan])}, 'x cannot contain NaN'), + ({'y': np.array([0.1, np.inf])}, 'y cannot contain NaN'), + ({'weights': np.array([1.0])}, 'Weights must have the same shape'), + ({'weights': np.array([1.0, np.nan])}, 'Weights cannot be NaN'), + ({'weights': np.array([1.0, 0.0])}, 'Weights must be strictly positive'), + ], + ) + def test_run_invalid_data(self, engine: DreamSampler, overrides, match) -> None: + """Shape mismatches and non-finite/non-positive data raise ValueError + before any sampling.""" + data = { + 'x': np.array([1.0, 2.0]), + 'y': np.array([0.1, 0.2]), + 'weights': np.array([1.0, 1.0]), + } + data.update(overrides) + with pytest.raises(ValueError, match=match): + engine.run(**data, samples=10, burn=0, thin=1) + + def test_run_rejects_none_weights(self, engine: DreamSampler) -> None: + """weights=None gets a clear ValueError instead of a shape error + from ``np.asarray(None)`` (CR-5).""" + with pytest.raises(ValueError, match='weights must not be None'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=None, + samples=10, + burn=0, + thin=1, + ) + + def test_run_basic(self, engine: DreamSampler, monkeypatch) -> None: + """Verify that run() returns a dict with expected keys.""" + mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) + + result = engine.run( + x=np.array([1.0, 2.0]), + y=np.array([0.1, 0.2]), + weights=np.array([1.0, 1.0]), + samples=100, + burn=20, + thin=2, + population=5, + ) + + assert isinstance(result, dict) + assert 'draws' in result + assert 'param_names' in result + assert 'internal_bumps_object' in result + assert 'logp' in result + mock_FitDriver.assert_called_once() + + def test_run_with_progress_callback(self, engine: DreamSampler, monkeypatch) -> None: + """Verify progress callback is wired up as a monitor.""" + mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) + progress_callback = MagicMock() + + result = engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=5, + thin=1, + progress_callback=progress_callback, + ) + + assert result is not None + call_kwargs = mock_FitDriver.call_args.kwargs + assert 'monitors' in call_kwargs + assert len(call_kwargs['monitors']) == 1 + assert isinstance(call_kwargs['monitors'][0], BumpsProgressMonitor) + + def test_run_aborted_by_user_raises_fit_error(self, engine: DreamSampler, monkeypatch) -> None: + """Verify that sampling abortion raises FitError.""" + self._setup_driver_mock(monkeypatch, fitter_state_value=TestDreamSamplerRun.ABORT) + + with pytest.raises(FitError, match='Sampling aborted by user'): + engine.run(x=np.array([1.0]), y=np.array([0.1]), weights=np.array([1.0])) + + def test_run_driver_exception_restores_parameters( + self, engine: DreamSampler, monkeypatch + ) -> None: + """Verify that a driver exception during sampling restores parameter values.""" + self._setup_driver_mock(monkeypatch, fit_side_effect=RuntimeError('driver failed')) + engine._restore_parameter_values = MagicMock() + + with pytest.raises(RuntimeError, match='driver failed'): + engine.run(x=np.array([1.0]), y=np.array([0.1]), weights=np.array([1.0])) + + engine._restore_parameter_values.assert_called_once() + + def test_run_population_param(self, engine: DreamSampler, monkeypatch) -> None: + """population kwarg is forwarded to DREAM as pop.""" + mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) + + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=0, + thin=1, + population=7, + ) + + call_kwargs = mock_FitDriver.call_args.kwargs + assert call_kwargs['pop'] == 7 + + def test_run_sampler_kwargs_forwarded(self, engine: DreamSampler, monkeypatch) -> None: + """sampler_kwargs entries are merged into the DREAM kwargs.""" + mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) + + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=0, + thin=1, + sampler_kwargs={'trim': False}, + ) + + assert mock_FitDriver.call_args.kwargs['trim'] is False + + def test_run_rejects_non_callable_callback(self, engine: DreamSampler) -> None: + with pytest.raises(ValueError, match='progress_callback must be callable'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=5, + thin=1, + progress_callback='not-callable', + ) + + # --- Resume-state tests ------------------------------------------------- + + def _make_resume_state_mock(self, *, nvar=2, npop=10, labels=None): + """Build a mock MCMCDraw for resume tests. + + BUMPS labels follow the pattern ``'p'`` (the + ``PARAMETER_PREFIX`` concatenated with the unique name), + e.g. ``'pFilm_thickness'``. + + :param nvar: Number of parameters. + :param npop: Population size. + :param labels: Parameter labels (defaults to ``['pa', 'pb']`` + which strip to ``['a', 'b']``). + """ + if labels is None: + labels = ['pa', 'pb'] + mock_state = MagicMock() + mock_state.Nvar = nvar + mock_state.Npop = npop + mock_state.labels = labels + mock_draw = MagicMock() + mock_draw.points = np.ones((20, nvar)) + mock_draw.logp = np.ones(20) + mock_state.draw.return_value = mock_draw + return mock_state + + def _make_problem_with_parameters(self, param_names): + """Build a mock FitProblem whose ``_parameters`` yields the given names.""" + params = [] + for name in param_names: + p = MagicMock() + p.name = 'p' + name + params.append(p) + mock_problem = MagicMock() + mock_problem._parameters = params + return mock_problem + + def test_run_resume_state(self, engine: DreamSampler, monkeypatch) -> None: + """Verify resume_state is forwarded to driver.fit().""" + mock_FitDriver, mock_driver = self._setup_driver_mock(monkeypatch) + resume_state = self._make_resume_state_mock() + self._set_problem(monkeypatch, self._make_problem_with_parameters(['a', 'b'])) + + result = engine.run( + x=np.array([1.0, 2.0]), + y=np.array([0.1, 0.2]), + weights=np.array([1.0, 1.0]), + samples=10, + burn=0, + thin=1, + resume_state=resume_state, + ) + + assert result is not None + # Verify a fit_state (defensive copy of resume_state) was passed to driver.fit() + call_kwargs = mock_driver.fit.call_args.kwargs + assert call_kwargs.get('fit_state') is not None + assert call_kwargs['fit_state'] is not resume_state + + def test_run_resume_param_mismatch_raises(self, engine: DreamSampler, monkeypatch) -> None: + """Parameter count mismatch raises ValueError before driver.fit().""" + self._set_problem(monkeypatch, self._make_problem_with_parameters(['a'])) + # resume_state has 2 params, model has 1 + resume_state = self._make_resume_state_mock(nvar=2) + + with pytest.raises(ValueError, match='resume_state has 2 parameters'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=0, + thin=1, + resume_state=resume_state, + ) + + def test_run_resume_param_name_mismatch_raises( + self, engine: DreamSampler, monkeypatch + ) -> None: + """Parameter name/order mismatch raises ValueError before driver.fit().""" + self._set_problem(monkeypatch, self._make_problem_with_parameters(['a', 'b'])) + # resume_state has labels ['px', 'py'] → stripped to ['x', 'y'] + # Current model has params ['pa', 'pb'] → stripped to ['a', 'b'] + # → mismatch + resume_state = self._make_resume_state_mock(nvar=2, labels=['px', 'py']) + + with pytest.raises(ValueError, match='Parameter names/order mismatch'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=0, + thin=1, + resume_state=resume_state, + ) + + def test_run_resume_population_mismatch_raises( + self, engine: DreamSampler, monkeypatch + ) -> None: + """Explicit population differing from state.Npop raises ValueError.""" + self._set_problem(monkeypatch, self._make_problem_with_parameters(['a', 'b'])) + resume_state = self._make_resume_state_mock(nvar=2, npop=10) + + with pytest.raises(ValueError, match='would produce'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=0, + thin=1, + population=3, # ceil(3*2)=6 ≠ 10 + resume_state=resume_state, + ) + + def test_run_resume_forces_burn_to_zero( + self, engine: DreamSampler, monkeypatch, caplog: 'pytest.LogCaptureFixture' + ) -> None: + """burn>0 with resume_state warns and is forced to 0.""" + mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) + self._set_problem(monkeypatch, self._make_problem_with_parameters(['a', 'b'])) + resume_state = self._make_resume_state_mock() + + with caplog.at_level(logging.WARNING, logger='easyscience.fitting.bumps'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=5, + thin=1, + resume_state=resume_state, + ) + + assert 'ignored on resume' in caplog.text + # burn must be forced to 0 in the kwargs passed to BUMPS + assert mock_FitDriver.call_args.kwargs['burn'] == 0 + + def test_run_resume_unlabeled_state_warns_and_uses_absolute_pop( + self, engine: DreamSampler, monkeypatch, caplog: 'pytest.LogCaptureFixture' + ) -> None: + """A state reloaded from disk carries default labels ('P0', ...), so + name validation is skipped with a warning, and the saved population is + reproduced as a negative pop (BUMPS' absolute-chain-count convention).""" + mock_FitDriver, _ = self._setup_driver_mock(monkeypatch) + self._set_problem(monkeypatch, self._make_problem_with_parameters(['a', 'b'])) + resume_state = self._make_resume_state_mock(nvar=2, npop=10, labels=['P0', 'P1']) + + with caplog.at_level(logging.WARNING, logger='easyscience.fitting.bumps'): + engine.run( + x=np.array([1.0]), + y=np.array([0.1]), + weights=np.array([1.0]), + samples=10, + burn=0, + thin=1, + resume_state=resume_state, + ) + + assert 'does not carry parameter names' in caplog.text + assert mock_FitDriver.call_args.kwargs['pop'] == -10 + + +class TestDreamSamplerProgressPayload: + """``DreamSampler._build_sample_progress_payload``.""" + + @pytest.fixture + def engine(self) -> DreamSampler: + return DreamSampler(obj='obj', fit_function='fit_function') + + def test_payload_structure_and_sampling_flag(self, engine: DreamSampler) -> None: + mock_problem = MagicMock() + mock_problem.chisq.side_effect = [25.0, 12.5] + mock_problem.labels.return_value = ['palpha'] + mock_problem.getp.return_value = np.array([1.0]) + + payload = engine._build_sample_progress_payload(mock_problem, 7, np.array([1.0]), 12.5) + + assert payload['iteration'] == 7 + assert payload['chi2'] == 25.0 + assert payload['reduced_chi2'] == 12.5 + assert payload['parameter_values'] == {'alpha': 1.0} + assert payload['sampling'] is True + assert payload['finished'] is False + assert payload['refresh_plots'] is False + # The nllf already computed by the sampler is reused — no model + # re-evaluation via setp. + mock_problem.chisq.assert_any_call(nllf=12.5, norm=False) + mock_problem.chisq.assert_any_call(nllf=12.5, norm=True) + mock_problem.setp.assert_not_called() + + def test_payload_keys(self, engine: DreamSampler) -> None: + """Same keys as the classical-fit payload, plus ``sampling``.""" + mock_problem = MagicMock() + mock_problem.chisq.side_effect = [10.0, 5.0] + mock_problem.labels.return_value = ['pa'] + mock_problem.getp.return_value = np.array([5.0]) + + payload = engine._build_sample_progress_payload(mock_problem, 1, np.array([5.0]), nllf=5.0) + + expected_keys = { + 'iteration', + 'chi2', + 'reduced_chi2', + 'parameter_values', + 'refresh_plots', + 'finished', + 'sampling', + } + assert set(payload.keys()) == expected_keys diff --git a/tests/unit/fitting/test_multi_fitter.py b/tests/unit/fitting/test_multi_fitter.py index ac897653..f47dfcb7 100644 --- a/tests/unit/fitting/test_multi_fitter.py +++ b/tests/unit/fitting/test_multi_fitter.py @@ -142,6 +142,33 @@ def test_handles_single_dataset(self): assert np.allclose(results[0].y_calc, [1.1, 2.1, 3.1]) +# =================================================================== +# MultiFitter._fit_function_wrapper +# =================================================================== + + +class TestFitFunctionWrapper: + def test_fit_function_restored_with_multiple_datasets(self): + """Wrapping must not leave ``fit_function`` pointing at the last + dataset's function (regression: sampling on a 2+ dataset MultiFitter + silently swapped the user-visible fit function).""" + fit_objects = [Line(1.0, 0.5), Line(2.0, 1.5)] + mf = MultiFitter(fit_objects, fit_objects) + original = mf.fit_function + assert original is fit_objects[0] + + x = [np.array([0.0, 1.0, 2.0]), np.array([0.0, 1.0])] + mf._dependent_dims = [(3,), (2,)] + wrapped = mf._fit_function_wrapper(x, flatten=True) + + assert mf.fit_function is original + + # Each wrapped section still evaluates its own dataset's function. + y = wrapped(np.zeros(5)) + expected = np.hstack([fit_objects[0](x[0]), fit_objects[1](x[1])]) + assert np.allclose(y, expected) + + # =================================================================== # MultiFitter._precompute_reshaping with weights=None # =================================================================== diff --git a/tests/unit/fitting/test_sampler.py b/tests/unit/fitting/test_sampler.py index 665c1d66..40f85ffc 100644 --- a/tests/unit/fitting/test_sampler.py +++ b/tests/unit/fitting/test_sampler.py @@ -50,8 +50,13 @@ def __init__(self, labels): self.labels = list(labels) -def _bumps_fitter_and_data(): - """Build a 2-parameter BUMPS MultiFitter over a small sine model.""" +def _fitter_and_data(): + """Build a 2-parameter MultiFitter over a small sine model. + + The fitter keeps its default (LMFit) minimizer: sampling no longer + requires switching to BUMPS, only an installed ``bumps`` package. + """ + pytest.importorskip('bumps') ref_sin = AbsSin(0.2, np.pi) sp = AbsSin(0.354, 3.05) sp.offset.fixed = False @@ -60,10 +65,6 @@ def _bumps_fitter_and_data(): y = ref_sin(x) weights = np.ones_like(x) f = MultiFitter([sp], [sp]) - try: - f.switch_minimizer('Bumps') - except AttributeError: - pytest.skip('BUMPS is not installed') return f, sp, x, y, weights @@ -209,31 +210,38 @@ def test_load_chain_rejects_bad_skip(self, tmp_path, skip): class TestSamplerErrorPaths: - def test_sample_requires_bumps(self): - """sample() must raise RuntimeError if the minimizer is not BUMPS — - and must not mutate the fitter (no needless minimizer rebuild).""" + def test_sample_requires_bumps_package(self, monkeypatch): + """sample() must raise RuntimeError when the bumps package is not + installed — regardless of the active minimizer — and must not touch + the fitter.""" sp = AbsSin(0.354, 3.05) f = MultiFitter([sp], [sp]) x, y, w = _xyw() sampler = Sampler(f, [x], [y], [w]) minimizer_before = f.minimizer - with pytest.raises(RuntimeError, match='Bayesian sampling requires a BUMPS minimizer'): + monkeypatch.setattr( + 'easyscience.fitting.available_minimizers.bumps_engine_available', False + ) + with pytest.raises(RuntimeError, match='requires the bumps package'): sampler.sample(samples=10, burn=5, thin=1) assert f.minimizer is minimizer_before - def test_fit_function_restored_on_error(self): - """fit_function must be restored even when the minimizer raises.""" - f, _, x, y, weights = _bumps_fitter_and_data() + def test_fitter_untouched_on_error(self): + """The fitter is never mutated by sampling, even when the engine + raises.""" + f, _, x, y, weights = _fitter_and_data() sampler = Sampler(f, [x], [y], [weights]) original_func = f.fit_function + minimizer_before = f.minimizer - # Invalid `samples` is rejected by the minimizer (single source of - # validation) *after* the fitter has been mutated for sampling. + # Invalid `samples` is rejected by the engine (single source of + # validation). with pytest.raises(ValueError, match='samples must be a positive integer'): sampler.sample(samples=-1, burn=5, thin=1) assert f.fit_function is original_func + assert f.minimizer is minimizer_before def test_extend_requires_existing_state(self): """extend() before sample()/load_state() raises RuntimeError.""" @@ -393,26 +401,13 @@ def test_fingerprint_without_weights(self): assert isinstance(sampler._fingerprint(), str) -class TestSamplingResultsLegacyDict: - def test_to_legacy_dict_maps_fields(self): - state = object() - results = SamplingResults( - draws=np.ones((2, 1)), param_names=['p'], logp=np.zeros(2), state=state - ) - legacy = results.to_legacy_dict() - assert legacy['internal_bumps_object'] is state - assert legacy['param_names'] == ['p'] - np.testing.assert_array_equal(legacy['draws'], results.draws) - np.testing.assert_array_equal(legacy['logp'], results.logp) - - class TestSamplerRunEngine: """The ``_run`` tail: results construction, storage, and kwarg merging, - with the minimizer's sampling entry point stubbed out.""" + with the ``DreamSampler`` engine stubbed out.""" def test_run_stores_results_and_exposes_properties(self, monkeypatch): - f, _, x, y, weights = _bumps_fitter_and_data() - from easyscience.fitting.minimizers.minimizer_bumps import Bumps + f, _, x, y, weights = _fitter_and_data() + from easyscience.fitting.samplers.sampler_dream import DreamSampler canned = { 'draws': np.arange(8.0).reshape(4, 2), @@ -422,14 +417,15 @@ def test_run_stores_results_and_exposes_properties(self, monkeypatch): } captured = {} - def fake_mcmc_sample(self, **kwargs): + def fake_run(self, **kwargs): captured.update(kwargs) return dict(canned) - monkeypatch.setattr(Bumps, 'mcmc_sample', fake_mcmc_sample) + monkeypatch.setattr(DreamSampler, 'run', fake_run) sampler = Sampler(f, [x], [y], [weights], sampler_kwargs={'trim': False}) original_func = f.fit_function + minimizer_before = f.minimizer results = sampler.sample(samples=100, burn=10, thin=2, sampler_kwargs={'init': 'lhs'}) assert isinstance(results, SamplingResults) @@ -443,8 +439,44 @@ def fake_mcmc_sample(self, **kwargs): assert captured['samples'] == 100 assert captured['burn'] == 10 assert captured['resume_state'] is None - # The fitter's fit function is restored after the run. + # The fitter is never mutated: a fresh engine gets the wrapped + # function directly, and the active (LMFit) minimizer stays put. assert f.fit_function is original_func + assert f.minimizer is minimizer_before + + def test_run_works_with_non_bumps_minimizer(self, monkeypatch): + """Sampling works with the default LMFit minimizer active — the + engine is constructed independently of the fitter's minimizer.""" + f, _, x, y, weights = _fitter_and_data() + from easyscience.fitting.samplers.sampler_dream import DreamSampler + + assert f.minimizer.package != 'bumps' # default is LMFit + + constructed = {} + original_init = DreamSampler.__init__ + + def spy_init(self, obj, fit_function): + constructed['obj'] = obj + constructed['fit_function'] = fit_function + original_init(self, obj, fit_function) + + canned = { + 'draws': np.zeros((2, 2)), + 'param_names': ['offset', 'phase'], + 'logp': np.zeros(2), + 'internal_bumps_object': object(), + } + monkeypatch.setattr(DreamSampler, '__init__', spy_init) + monkeypatch.setattr(DreamSampler, 'run', lambda self, **kwargs: dict(canned)) + + sampler = Sampler(f, [x], [y], [weights]) + results = sampler.sample(samples=10, burn=0, thin=1) + + assert results.param_names == ['offset', 'phase'] + # The engine is bound to the fitter's model object and a wrapped + # fit function, not to the minimizer. + assert constructed['obj'] is f.fit_object + assert callable(constructed['fit_function']) class TestSamplerExtendArithmetic: