diff --git a/docs/docs/tutorials/tutorial0_basics.ipynb b/docs/docs/tutorials/tutorial0_basics.ipynb index 10ec842f..7efb5705 100644 --- a/docs/docs/tutorials/tutorial0_basics.ipynb +++ b/docs/docs/tutorials/tutorial0_basics.ipynb @@ -416,7 +416,7 @@ " name='Straight line',\n", ")\n", "\n", - "binding = edyn.FitBinding(parameter_name='Gaussian area', model=fit_func)\n", + "binding = edyn.FitBinding(model=fit_func, targets='Gaussian area')\n", "\n", "parameter_analysis = edyn.ParameterAnalysis(\n", " parameters=analysis,\n", diff --git a/docs/docs/tutorials/tutorial0_more_advanced.ipynb b/docs/docs/tutorials/tutorial0_more_advanced.ipynb index 6869b259..da746957 100644 --- a/docs/docs/tutorials/tutorial0_more_advanced.ipynb +++ b/docs/docs/tutorials/tutorial0_more_advanced.ipynb @@ -339,11 +339,11 @@ " coefficients=[1.1, 0.2], x_unit='1/angstrom', y_unit='meV', name='DHO center fit'\n", ")\n", "\n", - "binding1 = edyn.FitBinding(parameter_name='Gaussian area', model=gauss_fit_func)\n", + "binding1 = edyn.FitBinding(model=gauss_fit_func, targets='Gaussian area')\n", "\n", - "binding2 = edyn.FitBinding(parameter_name='DHO area', model=dho_area_fit_func)\n", + "binding2 = edyn.FitBinding(model=dho_area_fit_func, targets='DHO area')\n", "\n", - "binding3 = edyn.FitBinding(parameter_name='DHO center', model=dho_center_fit_func)\n", + "binding3 = edyn.FitBinding(model=dho_center_fit_func, targets='DHO center')\n", "\n", "parameter_analysis = edyn.ParameterAnalysis(\n", " parameters=analysis,\n", diff --git a/docs/docs/tutorials/tutorial1_brownian.ipynb b/docs/docs/tutorials/tutorial1_brownian.ipynb index 269bb259..bb640325 100644 --- a/docs/docs/tutorials/tutorial1_brownian.ipynb +++ b/docs/docs/tutorials/tutorial1_brownian.ipynb @@ -451,7 +451,7 @@ "$$\n", "where $\\Gamma(Q) = D Q^2$ and $D$ is the diffusion coefficient. $S$ is an overall scale.\n", "\n", - "To fit the Brownian translational diffusion model to the data, we use the `ParameterAnalysis`. This time, we wish to fit the `Lorentzian`, and we wish to fit both its area (scale, $S$) and width to the diffusion model. We do this by creating a `FitBinding`, saying we want to fit both the area an width of the component called `Lorentzian`" + "To fit the Brownian translational diffusion model to the data, we use the `ParameterAnalysis`. This time, we wish to fit the `Lorentzian`, and we wish to fit both its area (scale, $S$) and width to the diffusion model. We do this by creating a `FitBinding`. Because the diffusion model's `lorentzian_name` matches the fitted component's name, the default targets automatically fit both the area and the width against the `Lorentzian area` and `Lorentzian width` parameters" ] }, { @@ -462,18 +462,17 @@ "outputs": [], "source": [ "brownian_diffusion_model = sm.BrownianTranslationalDiffusion(\n", - " name='Brownian Translational Diffusion', diffusion_coefficient=2.4e-9, scale=0.5\n", + " name='Brownian Translational Diffusion',\n", + " lorentzian_name='Lorentzian',\n", + " diffusion_coefficient=2.4e-9,\n", + " scale=0.5,\n", ")\n", "\n", - "binding = edyn.FitBinding(\n", - " parameter_name='Lorentzian',\n", - " model=brownian_diffusion_model,\n", - " modes=['area', 'width'],\n", - ")\n", + "binding = edyn.FitBinding(model=brownian_diffusion_model)\n", "\n", "parameter_analysis = edyn.ParameterAnalysis(\n", " parameters=diffusion_analysis,\n", - " bindings=[binding],\n", + " bindings=binding,\n", ")" ] }, @@ -495,16 +494,6 @@ "parameter_analysis.plot()" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "57b76d06", - "metadata": {}, - "outputs": [], - "source": [ - "parameter_analysis.plot(names=['Polynomial_c0'])" - ] - }, { "cell_type": "markdown", "id": "64babb01", @@ -691,7 +680,7 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": "default", "language": "python", "name": "python3" }, @@ -705,7 +694,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.14.5" + "version": "3.14.6" } }, "nbformat": 4, diff --git a/src/easydynamics/analysis/fit_binding.py b/src/easydynamics/analysis/fit_binding.py index da74e383..2b3cf7d3 100644 --- a/src/easydynamics/analysis/fit_binding.py +++ b/src/easydynamics/analysis/fit_binding.py @@ -3,6 +3,7 @@ from __future__ import annotations +from dataclasses import replace from typing import TYPE_CHECKING from easydynamics.base_classes.easydynamics_base import EasyDynamicsBase @@ -11,137 +12,112 @@ from easydynamics.sample_model.diffusion_model.diffusion_model_base import DiffusionModelBase if TYPE_CHECKING: - from collections.abc import Callable + from easydynamics.utils.fit_target import FitTarget class FitBinding(EasyDynamicsBase): """ - Contract between dataset, model, and fit function for ParameterAnalysis. This class - encapsulates the necessary information to bind a dataset key to a model and convert it into a - fit function callable. + Contract between dataset, model, and fit functions for ParameterAnalysis. A binding maps the + model's fittable predictions (its FitTargets) onto keys of the parameters Dataset they should + be fitted against. Examples -------- - **Basic usage with a ModelComponent** + **Fitting a component model to one parameter** - Bind a fit parameter name to any ``ModelComponent`` or ``ComponentCollection``: + Component models (e.g. a Polynomial) have a single prediction — their ``evaluate`` — so + ``targets`` is simply the dataset key to fit against. The model's x_unit/y_unit declare the + units its evaluate expects: here x is the dataset's Q coordinate and y the fitted parameter, so + construct the model with matching units (or pass ``x_unit=None`` / ``y_unit=None`` to fit raw + values): ```python import easydynamics as edyn import easydynamics.sample_model as sm - fit_func = sm.Polynomial(coefficients=[3.7, -0.5], display_name='Straight line') - binding = edyn.FitBinding(parameter_name='Gaussian area', model=fit_func) + fit_func = sm.Polynomial( + coefficients=[3.7, -0.5], + x_unit='1/angstrom', + y_unit='meV', + display_name='Straight line', + ) + binding = edyn.FitBinding(model=fit_func, targets='Gaussian area') ``` - **Usage with a DiffusionModelBase and specific modes** + **Fitting a diffusion model with default dataset keys** - For diffusion models, use ``modes`` to select which parameter arrays are fitted: + Diffusion models declare their predictions (``'area'``, ``'width'``, and for DeltaLorentz also + ``'delta_area'``). With ``targets=None`` all predictions are fitted against default dataset + keys derived from the model's component names: ```python brownian = sm.BrownianTranslationalDiffusion( - display_name='Brownian Translational Diffusion', diffusion_coefficient=2.4e-9, scale=0.5, + lorentzian_name='Lorentzian', ) + binding = edyn.FitBinding(model=brownian) # fits 'Lorentzian area' and 'Lorentzian width' + ``` + + **Selecting predictions or mapping them to custom dataset keys** + + Pass a list of prediction names, or a dict mapping prediction names to dataset keys: + ```python + binding = edyn.FitBinding(model=brownian, targets=['width']) + + delta_lorentz = sm.DeltaLorentz(A_0=0.5, lorentzian_width=0.1) binding = edyn.FitBinding( - parameter_name='Lorentzian', - model=brownian, - modes=['area', 'width'], + model=delta_lorentz, + targets={ + 'width': 'Lorentzian width', + 'area': 'Lorentzian area', + 'delta_area': 'Elastic area', + }, ) ``` """ def __init__( self, - parameter_name: str, model: ModelComponent | ComponentCollection | DiffusionModelBase, - modes: str | list[str] | None = None, + targets: str | list[str] | dict[str, str] | None = None, display_name: str | None = None, unique_name: str | None = None, ) -> None: """ Initialize a FitBinding. + Validation raises ``TypeError`` if model or targets have an invalid type, and + ``ValueError`` if targets names a prediction the model does not declare. + Parameters ---------- - parameter_name : str - The name of the parameter to fit. This should correspond to a key in the dataset. model : ModelComponent | ComponentCollection | DiffusionModelBase The model to fit. This can be a single ModelComponent, a ComponentCollection, or a DiffusionModelBase. - modes : str | list[str] | None, default=None - The modes to fit for diffusion models. This can be a single string, a list of strings, - or None (which defaults to ["area", "width"]). Only applicable if the model is a - DiffusionModelBase. Default is None. + targets : str | list[str] | dict[str, str] | None, default=None + Which predictions of the model to fit, and against which dataset keys. For component + models this must be a string: the dataset key to fit the model's ``evaluate`` against. + For diffusion models: None fits all predictions against their default dataset keys; a + string or list of strings selects predictions by name (default keys); a dict maps + prediction names to custom dataset keys. display_name : str | None, default=None An optional display name for the FitBinding. If None, the unique_name will be used. Default is None. unique_name : str | None, default=None An optional unique name for the FitBinding. If None, a unique name will be generated. Default is None. - - Raises - ------ - TypeError - If parameter_name is not a string, if model is not a ModelComponent, - ComponentCollection or DiffusionModelBase, or if modes is not a string, list of - strings, or None. """ super().__init__(display_name=display_name, unique_name=unique_name) - if not isinstance(parameter_name, str): - raise TypeError('parameter_name must be a string') - - if not isinstance(model, (ModelComponent, ComponentCollection, DiffusionModelBase)): - raise TypeError( - 'model must be a ModelComponent, ComponentCollection, or DiffusionModelBase' - ) - - if modes is not None and not isinstance(modes, (str, list)): - raise TypeError('modes must be a string, list of strings, or None') - - if isinstance(modes, list) and not all(isinstance(mode, str) for mode in modes): - raise TypeError('All modes in the list must be strings') - - self._parameter_name = parameter_name + self._validate_model(model) + self._normalize_targets(model, targets) self._model = model - self._modes = modes + self._targets = targets # ------------------------------------------------------------------ # Properties # ------------------------------------------------------------------ - @property - def parameter_name(self) -> str: - """ - The name of the parameter to fit. This should correspond to a key in the dataset. - - Returns - ------- - str - The name of the parameter to fit. - """ - return self._parameter_name - - @parameter_name.setter - def parameter_name(self, value: str) -> None: - """ - Set the name of the parameter to fit. - - Parameters - ---------- - value : str - The new name of the parameter to fit. - - Raises - ------ - TypeError - If the value is not a string. - """ - if not isinstance(value, str): - raise TypeError('parameter_name must be a string') - self._parameter_name = value - @property def model(self) -> ModelComponent | ComponentCollection | DiffusionModelBase: """ @@ -160,162 +136,164 @@ def model(self, value: ModelComponent | ComponentCollection | DiffusionModelBase """ Set the model to fit. + Validation raises ``TypeError`` if the value has an invalid type, and ``ValueError`` if the + current targets name a prediction the new model does not declare. + Parameters ---------- value : ModelComponent | ComponentCollection | DiffusionModelBase The new model to fit. - - Raises - ------ - TypeError - If the value is not a ModelComponent, ComponentCollection, or DiffusionModelBase. """ - if not isinstance(value, (ModelComponent, ComponentCollection, DiffusionModelBase)): - raise TypeError( - 'model must be a ModelComponent, ComponentCollection, or DiffusionModelBase.' - ) + self._validate_model(value) + self._normalize_targets(value, self._targets) self._model = value @property - def modes(self) -> str | list[str] | None: + def targets(self) -> str | list[str] | dict[str, str] | None: """ - The modes to fit for diffusion models. This can be a single string, a list of strings, or - None (which defaults to ["area", "width"]). + Which predictions of the model to fit, and against which dataset keys. Returns ------- - str | list[str] | None - The modes to fit for diffusion models. + str | list[str] | dict[str, str] | None + The targets specification (see ``__init__``). """ - return self._modes + return self._targets - @modes.setter - def modes(self, value: str | list[str] | None) -> None: + @targets.setter + def targets(self, value: str | list[str] | dict[str, str] | None) -> None: """ - Set the modes to fit for diffusion models. + Set which predictions of the model to fit, and against which dataset keys. + + Validation raises ``TypeError`` if the value has an invalid type for the current model, and + ``ValueError`` if it names a prediction the model does not declare. Parameters ---------- - value : str | list[str] | None - The new modes to fit for diffusion models. - - Raises - ------ - TypeError - If the value is not a string, list of strings, or None. + value : str | list[str] | dict[str, str] | None + The new targets specification (see ``__init__``). """ - if value is not None and not isinstance(value, (str, list)): - raise TypeError('modes must be a string, list of strings, or None') - - if isinstance(value, str): - value = [value] - if isinstance(value, list) and not all(isinstance(mode, str) for mode in value): - raise TypeError('All modes in the list must be strings') - self._modes = value + self._normalize_targets(self._model, value) + self._targets = value # ------------------------------------------------------------------ # Other methods # ------------------------------------------------------------------ - def build_callables(self) -> list[Callable]: - """ - Build the callables for fitting based on the model and modes. - - Returns - ------- - list[Callable] - A list of callables for fitting. - """ - modes = self._get_modes() - - if isinstance(self.model, DiffusionModelBase): - return [self._build_diffusion_callable(mode) for mode in modes] - - return [lambda x, **_: self.model.evaluate(x)] - - def get_model_names(self) -> list[str]: - """ - Get the names of the models based on the current modes. - - Returns - ------- - list[str] - A list of model names. + def get_targets(self) -> list[FitTarget]: """ - modes = self._get_modes() - - if isinstance(self.model, DiffusionModelBase): - return [f'{self.model.display_name} {mode}' for mode in modes] + Get the FitTargets this binding fits, with dataset keys resolved. - return [self.model.display_name] - - def get_parameter_names(self) -> list[str]: - """ - Get the names of the parameters based on the current modes. + Targets are built from the model at call time, so their units and default dataset keys + reflect the model's current state. Returns ------- - list[str] - A list of parameter names. + list[FitTarget] + The resolved fit targets. """ - modes = self._get_modes() - - if isinstance(self.model, DiffusionModelBase): - # This needs to be generalised. - # TODO: Generalise this for different diffusion models and modes. # noqa TD002 TD003 - if 'delta' in modes: - return [f'{self.parameter_name} area' for mode in modes] - - return [f'{self.parameter_name} {mode}' for mode in modes] - - return [self.parameter_name] + available = {target.name: target for target in self.model.get_fit_targets()} + normalized = self._normalize_targets(self.model, self._targets) + return [ + available[name] + if dataset_key is None + else replace(available[name], dataset_key=dataset_key) + for name, dataset_key in normalized.items() + ] # ------------------------------------------------------------------ # Private methods # ------------------------------------------------------------------ - def _build_diffusion_callable(self, mode: str) -> Callable: + @staticmethod + def _validate_model( + model: ModelComponent | ComponentCollection | DiffusionModelBase, + ) -> None: """ - Build a callable for a specific diffusion mode. + Validate the model type. Parameters ---------- - mode : str - The diffusion mode ("area" or "width"). - - Returns - ------- - Callable - A callable for the specified diffusion mode. + model : ModelComponent | ComponentCollection | DiffusionModelBase + The model to validate. Raises ------ - ValueError - If the mode is unknown. + TypeError + If model is not a ModelComponent, ComponentCollection, or DiffusionModelBase. """ - model = self.model - - if mode == 'area': - return lambda x, **_: model.calculate_QISF(x) * model.scale.value - - if mode == 'width': - return lambda x, **_: model.calculate_width(x) + if not isinstance(model, (ModelComponent, ComponentCollection, DiffusionModelBase)): + raise TypeError( + 'model must be a ModelComponent, ComponentCollection, or DiffusionModelBase' + ) - if mode == 'delta': - return lambda x, **_: model.calculate_EISF(x) * model.scale.value + @staticmethod + def _normalize_targets( + model: ModelComponent | ComponentCollection | DiffusionModelBase, + targets: str | list[str] | dict[str, str] | None, + ) -> dict[str, str | None]: + """ + Validate a targets specification and normalize it to a prediction-name mapping. - raise ValueError(f'Unknown diffusion mode: {mode}') + This is the one place that knows the two spec dialects: component models take the dataset + key as a plain string for their single ``'value'`` prediction, while diffusion models + select predictions by name (None selects all) with optional dataset-key overrides. - def _get_modes(self) -> list[str]: - """ - Get the modes to fit for diffusion models, defaulting to ["area", "width"] if not set. + Parameters + ---------- + model : ModelComponent | ComponentCollection | DiffusionModelBase + The model the targets apply to. + targets : str | list[str] | dict[str, str] | None + The targets specification to validate and normalize. Returns ------- - list[str] - The modes to fit for diffusion models. + dict[str, str | None] + Mapping of prediction name to dataset-key override (None means the prediction's default + dataset key is used). + + Raises + ------ + TypeError + If targets has an invalid type for the given model. + ValueError + If targets names a prediction the model does not declare. """ - return ['area', 'width'] if self.modes is None else self.modes + if not isinstance(model, DiffusionModelBase): + if not isinstance(targets, str): + raise TypeError( + 'For component models, targets must be the dataset key (a string) to fit ' + "the model's evaluate against." + ) + return {'value': targets} + + available = [target.name for target in model.get_fit_targets()] + if targets is None: + return dict.fromkeys(available) + if isinstance(targets, str): + requested = [targets] + elif isinstance(targets, list): + requested = targets + elif isinstance(targets, dict): + requested = list(targets.keys()) + if not all(isinstance(key, str) for key in targets.values()): + raise TypeError('targets dict values must be dataset keys (strings)') + else: + raise TypeError( + 'targets must be None, a prediction name, a list of prediction names, ' + 'or a dict mapping prediction names to dataset keys' + ) + if not all(isinstance(name, str) for name in requested): + raise TypeError('prediction names in targets must be strings') + + unknown = sorted(set(requested) - set(available)) + if unknown: + raise ValueError( + f'Unknown prediction(s) {", ".join(unknown)} for ' + f'{model.__class__.__name__}. Available predictions: ' + f'{", ".join(available)}.' + ) + return dict(targets) if isinstance(targets, dict) else dict.fromkeys(requested) # ------------------------------------------------------------------ # dunder methods @@ -331,10 +309,8 @@ def __repr__(self) -> str: A string representation of the FitBinding. """ return ( - f'{self.__class__.__name__}(' - f'parameter_name={self.parameter_name!r},\n' - f' model={self.model.display_name!r},\n' - f' modes={self.modes},\n' - f' display_name={self.display_name!r},\n' - f' unique_name={self.unique_name!r})' + f'{self.__class__.__name__}(\n' + f' model={self.model.display_name},\n' + f' targets={self.targets},\n' + f')' ) diff --git a/src/easydynamics/analysis/parameter_analysis.py b/src/easydynamics/analysis/parameter_analysis.py index 142f0b1a..7e24108e 100644 --- a/src/easydynamics/analysis/parameter_analysis.py +++ b/src/easydynamics/analysis/parameter_analysis.py @@ -15,7 +15,9 @@ from easydynamics.analysis.analysis import Analysis from easydynamics.analysis.fit_binding import FitBinding from easydynamics.base_classes.easydynamics_modelbase import EasyDynamicsModelBase +from easydynamics.utils.fit_target import FitTarget from easydynamics.utils.utils import _in_notebook +from easydynamics.utils.utils import convert_value_unit class ParameterAnalysis(EasyDynamicsModelBase): @@ -25,15 +27,15 @@ class ParameterAnalysis(EasyDynamicsModelBase): Can be used to fit parameters to ModelComponents, ComponentCollections, or DiffusionModelBase objects, and to plot the parameters and fit results. The parameters to be analyzed can be provided as a sc.Dataset or directly as an Analysis object. Multiple parameters can be fitted - simultaneously, and the fit functions can be customized for each parameter. For diffusion - models, the area and width can be fitted separately (or not at all) by specifying fit settings. + simultaneously, and each binding maps its model's predictions onto the dataset keys they are + fitted against (for diffusion models e.g. 'area', 'width', or 'delta_area'). Examples -------- **Fitting Lorentzian widths to a diffusion model** - After a full Analysis fit, pass the Analysis directly and bind each parameter to a fit function - using a ``FitBinding``: + After a full Analysis fit, pass the Analysis directly and bind the model's predictions to + dataset keys using a ``FitBinding``: ```python import easydynamics as edyn import easydynamics.sample_model as sm @@ -41,9 +43,8 @@ class ParameterAnalysis(EasyDynamicsModelBase): # analysis is an edyn.Analysis object with previously fitted parameters diffusion_model = sm.BrownianTranslationalDiffusion(diffusion_coefficient=2.4e-9, scale=0.5) binding = edyn.FitBinding( - parameter_name='Lorentzian width', model=diffusion_model, - modes=['width'], + targets={'width': 'Lorentzian width'}, ) param_analysis = edyn.ParameterAnalysis( @@ -56,10 +57,13 @@ class ParameterAnalysis(EasyDynamicsModelBase): **Fitting multiple parameters with separate bindings** + Component models declare the units their evaluate expects: here the Polynomial's x is the + dataset's Q coordinate and its y is the fitted parameter, so construct it with matching units + (or pass ``x_unit=None`` / ``y_unit=None`` to fit raw values): ```python area_binding = edyn.FitBinding( - parameter_name='Lorentzian area', - model=sm.Polynomial(coefficients=[0.5, 0.0]), + model=sm.Polynomial(coefficients=[0.5, 0.0], x_unit='1/angstrom', y_unit='meV'), + targets='Lorentzian area', ) param_analysis = edyn.ParameterAnalysis( parameters=analysis, @@ -183,23 +187,24 @@ def fit(self) -> FitResults: funcs, models = [], [] for binding in self.bindings: - param_names = binding.get_parameter_names() - callables = binding.build_callables() - - for pname, func in zip(param_names, callables, strict=True): - if pname not in self.parameters: + for target in binding.get_targets(): + if target.dataset_key not in self.parameters: raise ValueError( - f"Parameter '{pname}' from binding '{binding.unique_name}' " - f'not found in parameters Dataset.' + f"Parameter '{target.dataset_key}' from binding " + f"'{binding.unique_name}' not found in parameters Dataset." ) - x, y, weight = self._get_xyweight_from_dataset(pname) + x, y, weight = self._get_xyweight_from_dataset(target.dataset_key) + x_factor, y_factor = self._get_unit_conversions(target) + x = x * x_factor + y = y * y_factor + weight = weight / y_factor xs.append(x) ys.append(y) ws.append(weight) - funcs.append(func) + funcs.append(target.function) models.append(binding.model) mf = MultiFitter( @@ -258,7 +263,7 @@ def plot( names = list(self.parameters.keys()) else: for b in self.bindings: - names.extend(b.get_parameter_names()) + names.extend(target.dataset_key for target in b.get_targets()) names = self._normalize_names(names) @@ -282,14 +287,13 @@ def plot( data_arrays = {} model_arrays = {} - # map parameter names to model names + # map dataset keys to model labels param_to_model = {} if self.bindings is not None: for b in self.bindings: - param_names = b.get_parameter_names() - model_names = b.get_model_names() - - param_to_model.update(dict(zip(param_names, model_names, strict=True))) + param_to_model.update({ + target.dataset_key: target.label for target in b.get_targets() + }) for pname in names: data_arrays[pname] = self.parameters[pname] @@ -354,22 +358,20 @@ def calculate_model_dataset(self, bindings: list[FitBinding]) -> sc.Dataset: arrays = {} for b in bindings: - param_names = b.get_parameter_names() - model_names = b.get_model_names() - callables = b.build_callables() - - for pname, mname, func in zip(param_names, model_names, callables, strict=True): - if pname not in self.parameters: + for target in b.get_targets(): + if target.dataset_key not in self.parameters: raise ValueError( - f"Parameter '{pname}' from binding '{b.unique_name}' " + f"Parameter '{target.dataset_key}' from binding '{b.unique_name}' " f'not found in parameters Dataset.' ) - da = self.parameters[pname] + da = self.parameters[target.dataset_key] x = da.coords['Q'] + x_factor, y_factor = self._get_unit_conversions(target) - y_model = func(x.values) + y_model = target.function(x.values * x_factor) + y_model = y_model / y_factor - arrays[mname] = sc.DataArray( + arrays[target.label] = sc.DataArray( data=sc.array(dims=['Q'], values=y_model, unit=da.unit), coords={'Q': x}, ) @@ -521,6 +523,55 @@ def _normalize_names(self, names: str | list[str] | None) -> list[str] | None: # Private methods ############# + def _get_unit_conversions(self, target: FitTarget) -> tuple[float, float]: + """ + Return (x_factor, y_factor) to convert dataset values into the target's declared units. + + x_factor converts Q coordinate values from their stored unit to target.x_unit (e.g. + 1/angstrom for diffusion-model predictions). y_factor converts parameter values from their + stored unit to target.y_unit. A factor is 1.0 when the target declares the corresponding + unit as None (its function takes raw values). ``sc.UnitError`` is raised when x or y units + are physically incompatible (e.g. meV vs 1/angstrom). + + Parameters + ---------- + target : FitTarget + The fit target whose unit contract defines the target units. + + Returns + ------- + tuple[float, float] + ``(x_factor, y_factor)`` scale factors to apply before/after model evaluation. + """ + da = self.parameters[target.dataset_key] + + def factor(from_unit: str, to_unit: str | None, error_message: str) -> float: + if to_unit is None: + return 1.0 + try: + return convert_value_unit(1.0, from_unit, str(to_unit)) + except Exception as e: + raise sc.UnitError(error_message) from e + + q_unit = str(da.coords['Q'].unit) + x_factor = factor( + q_unit, + target.x_unit, + f"Q coordinate unit '{q_unit}' is incompatible with " + f"the x_unit '{target.x_unit}' of fit target '{target.label}' " + f"for parameter '{target.dataset_key}'.", + ) + + param_unit = str(da.unit) + y_factor = factor( + param_unit, + target.y_unit, + f"Parameter '{target.dataset_key}' unit '{param_unit}' is incompatible " + f"with the y_unit '{target.y_unit}' of fit target '{target.label}'.", + ) + + return x_factor, y_factor + def _get_xyweight_from_dataset( self, parameter_name: str ) -> tuple[np.ndarray, np.ndarray, np.ndarray]: @@ -549,21 +600,29 @@ def _get_xyweight_from_dataset( raise ValueError(f"Parameter name '{parameter_name}' not found in parameters Dataset.") variances = self._parameters[parameter_name].variances + values = self._parameters[parameter_name].values + q_values = self._parameters[parameter_name].coords['Q'].values + if variances is None: - weight = np.ones_like(self._parameters[parameter_name].values) - elif np.any(~np.isfinite(variances)) or np.any(variances <= 0): + return q_values, values, np.ones_like(values) + + # NaN variances arise when a parameter is absent for a given Q (parameters_to_dataset + # fills np.nan for missing parameters). Filter those rows silently; other non-finite or + # non-positive variances are errors. + nan_mask = np.isnan(variances) + if np.any(~nan_mask & (~np.isfinite(variances) | (variances <= 0))): raise ValueError( f"Non-finite variances found for parameter '{parameter_name}', " f'cannot compute weights.' ) - else: - weight = 1 / np.sqrt(variances) + valid_mask = ~nan_mask + if not np.any(valid_mask): + raise ValueError( + f"No finite positive variances found for parameter '{parameter_name}', " + f'cannot compute weights.' + ) - return ( - self._parameters[parameter_name].coords['Q'].values, - self._parameters[parameter_name].values, - weight, - ) + return q_values[valid_mask], values[valid_mask], 1 / np.sqrt(variances[valid_mask]) ############# # Dunder methods @@ -587,9 +646,8 @@ def __repr__(self) -> str: binding_info = [ { - 'parameter': b.parameter_name, 'model': b.model.display_name, - 'modes': b.modes, + 'targets': b.targets, } for b in self._bindings ] diff --git a/tests/unit/easydynamics/analysis/test_fit_binding.py b/tests/unit/easydynamics/analysis/test_fit_binding.py index 191d1ba1..901ab9e5 100644 --- a/tests/unit/easydynamics/analysis/test_fit_binding.py +++ b/tests/unit/easydynamics/analysis/test_fit_binding.py @@ -2,288 +2,272 @@ # SPDX-License-Identifier: BSD-3-Clause -from unittest.mock import Mock - +import numpy as np import pytest from easydynamics.analysis.fit_binding import FitBinding from easydynamics.sample_model.components.gaussian import Gaussian +from easydynamics.sample_model.components.polynomial import Polynomial from easydynamics.sample_model.diffusion_model.brownian_translational_diffusion import ( BrownianTranslationalDiffusion, ) +from easydynamics.sample_model.diffusion_model.delta_lorentz import DeltaLorentz +from easydynamics.utils.fit_target import FitTarget class TestFitBinding: @pytest.fixture - def fit_binding(self): - model = Gaussian() - return FitBinding(parameter_name='parameter1', model=model) + def component_binding(self): + model = Gaussian(display_name='GaussianModel') + return FitBinding(model=model, targets='parameter1') @pytest.fixture - def diffusion_fit_binding(self): - model = BrownianTranslationalDiffusion() - return FitBinding(parameter_name='parameter3', model=model) - - def test_initialization(self, fit_binding): - # WHEN THEN EXPECT - assert isinstance(fit_binding, FitBinding) - assert fit_binding.parameter_name == 'parameter1' - assert isinstance(fit_binding.model, Gaussian) - assert fit_binding.modes is None - - @pytest.mark.parametrize( - 'parameter_name, model, modes, error_msg', - [ - # parameter_name errors - (123, Gaussian(), None, 'parameter_name must be a string'), - (None, Gaussian(), None, 'parameter_name must be a string'), - # model errors - ( - 'param', - 123, - None, - 'model must be a ModelComponent, ComponentCollection, or DiffusionModelBase', - ), - ( - 'param', - 'not_a_model', - None, - 'model must be a ModelComponent, ComponentCollection, or DiffusionModelBase', - ), - # modes type errors - ( - 'param', - Gaussian(), - 123, - 'modes must be a string, list of strings, or None', - ), - ( - 'param', - Gaussian(), - {'mode': 'area'}, - 'modes must be a string, list of strings, or None', - ), - # modes list content errors - ( - 'param', - Gaussian(), - ['area', 123], - 'All modes in the list must be strings', - ), - ('param', Gaussian(), [None], 'All modes in the list must be strings'), - ], - ) - def test_fitbinding_init_errors(self, parameter_name, model, modes, error_msg): - with pytest.raises(TypeError, match=error_msg): - FitBinding( - parameter_name=parameter_name, - model=model, - modes=modes, - ) + def diffusion_binding(self): + model = BrownianTranslationalDiffusion(lorentzian_name='Lorentzian') + return FitBinding(model=model) # ------------------------------------------------------------------ - # Properties + # Initialization and validation # ------------------------------------------------------------------ - def test_parameter_name_setter(self, fit_binding): - # WHEN - fit_binding.parameter_name = 'new_parameter' + def test_initialization(self, component_binding): + # WHEN THEN EXPECT + assert isinstance(component_binding, FitBinding) + assert component_binding.targets == 'parameter1' + assert isinstance(component_binding.model, Gaussian) - # THEN EXPECT - assert fit_binding.parameter_name == 'new_parameter' + def test_initialization_invalid_model_raises(self): + # WHEN THEN EXPECT + with pytest.raises(TypeError, match='model must be a ModelComponent'): + FitBinding(model='not_a_model', targets='parameter1') - def test_parameter_name_setter_errors(self, fit_binding): - with pytest.raises(TypeError, match='parameter_name must be a string'): - fit_binding.parameter_name = 123 + @pytest.mark.parametrize( + 'targets', + [None, ['parameter1'], {'value': 'parameter1'}, 123], + ids=['none', 'list', 'dict', 'int'], + ) + def test_component_model_requires_string_targets(self, targets): + # WHEN THEN EXPECT: component models have a single prediction, so targets must be + # the dataset key + with pytest.raises(TypeError, match='targets must be the dataset key'): + FitBinding(model=Gaussian(), targets=targets) - def test_model_setter(self, fit_binding): - # WHEN + def test_diffusion_unknown_prediction_raises(self): + # WHEN THEN EXPECT model = BrownianTranslationalDiffusion() + with pytest.raises(ValueError, match=r'Unknown prediction.*Available predictions'): + FitBinding(model=model, targets=['nonsense']) - # THEN - fit_binding.model = model - - # EXPECT - assert fit_binding.model is model + def test_diffusion_unknown_prediction_in_dict_raises(self): + # WHEN THEN EXPECT + model = BrownianTranslationalDiffusion() + with pytest.raises(ValueError, match='Unknown prediction'): + FitBinding(model=model, targets={'delta_area': 'Elastic area'}) - def test_model_setter_errors(self, fit_binding): - with pytest.raises( - TypeError, - match='model must be a ModelComponent, ComponentCollection, or DiffusionModelBase', - ): - fit_binding.model = 'not_a_model' + def test_diffusion_invalid_targets_type_raises(self): + # WHEN THEN EXPECT + model = BrownianTranslationalDiffusion() + with pytest.raises(TypeError, match='targets must be None'): + FitBinding(model=model, targets=123) - def test_modes_setter(self, fit_binding): - # WHEN - fit_binding.modes = 'area' + def test_diffusion_non_string_dataset_key_raises(self): + # WHEN THEN EXPECT + model = BrownianTranslationalDiffusion() + with pytest.raises(TypeError, match='dataset keys'): + FitBinding(model=model, targets={'width': 123}) - # THEN EXPECT - assert fit_binding.modes == ['area'] + # ------------------------------------------------------------------ + # Properties + # ------------------------------------------------------------------ - # WHEN - fit_binding.modes = ['area', 'width'] + def test_model_setter_revalidates_targets(self): + # WHEN: a binding using DeltaLorentz's delta_area prediction + binding = FitBinding(model=DeltaLorentz(), targets=['delta_area']) - # THEN EXPECT - assert fit_binding.modes == ['area', 'width'] + # THEN EXPECT: switching to a model without that prediction fails + with pytest.raises(ValueError, match='Unknown prediction'): + binding.model = BrownianTranslationalDiffusion() - def test_modes_setter_errors(self, fit_binding): - with pytest.raises(TypeError, match='modes must be a string, list of strings, or None'): - fit_binding.modes = 123 + def test_model_setter_invalid_type_raises(self, component_binding): + # WHEN THEN EXPECT + with pytest.raises(TypeError, match='model must be a ModelComponent'): + component_binding.model = 'not_a_model' - with pytest.raises(TypeError, match='modes must be a string, list of strings, or None'): - fit_binding.modes = {'mode': 'area'} + def test_targets_setter(self, diffusion_binding): + # WHEN THEN + diffusion_binding.targets = ['width'] - with pytest.raises(TypeError, match='All modes in the list must be strings'): - fit_binding.modes = ['area', 123] + # EXPECT + assert [t.name for t in diffusion_binding.get_targets()] == ['width'] - with pytest.raises(TypeError, match='All modes in the list must be strings'): - fit_binding.modes = [None] + def test_targets_setter_invalid_raises(self, diffusion_binding): + # WHEN THEN EXPECT + with pytest.raises(ValueError, match='Unknown prediction'): + diffusion_binding.targets = ['nonsense'] # ------------------------------------------------------------------ - # Other methods + # get_targets # ------------------------------------------------------------------ - def test_build_callables_component(self, fit_binding): + def test_component_target(self, component_binding): # WHEN - mock_model = Mock() - mock_model.evaluate.return_value = 1.0 - fit_binding._model = mock_model - - # THEN - callables = fit_binding.build_callables() + targets = component_binding.get_targets() + + # THEN EXPECT: a single target evaluating the component, with its units + assert len(targets) == 1 + target = targets[0] + assert isinstance(target, FitTarget) + assert target.name == 'value' + assert target.dataset_key == 'parameter1' + assert target.label == 'GaussianModel' + assert target.x_unit == component_binding.model.x_unit + assert target.y_unit == component_binding.model.y_unit + + def test_component_target_function_evaluates_model(self, component_binding): + # WHEN + target = component_binding.get_targets()[0] + x = np.linspace(-1, 1, 5) - # EXPECT - assert len(callables) == 1 - assert callable(callables[0]) - assert callables[0](0) == pytest.approx(1.0) - mock_model.evaluate.assert_called_once_with(0) + # THEN EXPECT + np.testing.assert_allclose(target.function(x), component_binding.model.evaluate(x)) - def test_build_callables_diffusion(self, diffusion_fit_binding): + def test_diffusion_default_targets(self, diffusion_binding): # WHEN - mock_model = Mock(spec=BrownianTranslationalDiffusion) - mock_model.calculate_QISF.return_value = 2.0 - mock_model.scale.value = 3.0 - mock_model.calculate_width.return_value = 0.5 - diffusion_fit_binding._model = mock_model + targets = diffusion_binding.get_targets() - # THEN - callables = diffusion_fit_binding.build_callables() + # THEN EXPECT: all predictions with default dataset keys from the Lorentzian name + assert [t.name for t in targets] == ['area', 'width'] + assert [t.dataset_key for t in targets] == ['Lorentzian area', 'Lorentzian width'] - # EXPECT - assert len(callables) == 2 - assert callable(callables[0]) - assert callable(callables[1]) - assert callables[0](0) == pytest.approx(6.0) # 2.0 * 3.0 - assert callables[1](0) == pytest.approx(0.5) - mock_model.calculate_QISF.assert_called_once_with(0) - mock_model.calculate_width.assert_called_once_with(0) - - def test_build_callables_diffusion_with_modes(self, diffusion_fit_binding): + def test_diffusion_string_target(self): # WHEN - diffusion_fit_binding.modes = 'area' - mock_model = Mock(spec=BrownianTranslationalDiffusion) - mock_model.calculate_QISF.return_value = 2.0 - mock_model.scale.value = 3.0 - diffusion_fit_binding._model = mock_model + model = BrownianTranslationalDiffusion(lorentzian_name='QENS') + binding = FitBinding(model=model, targets='width') # THEN - callables = diffusion_fit_binding.build_callables() + targets = binding.get_targets() # EXPECT - assert len(callables) == 1 - assert callable(callables[0]) - assert callables[0](0) == pytest.approx(6.0) # 2.0 * 3.0 - mock_model.calculate_QISF.assert_called_once_with(0) + assert len(targets) == 1 + assert targets[0].name == 'width' + assert targets[0].dataset_key == 'QENS width' - def test_get_model_names(self, fit_binding): - # WHEN THEN - model_names = fit_binding.get_model_names() - - # EXPECT - assert model_names == ['Gaussian'] + def test_diffusion_list_targets(self): + # WHEN + model = BrownianTranslationalDiffusion(lorentzian_name='QENS') + binding = FitBinding(model=model, targets=['width', 'area']) + + # THEN EXPECT: order follows the user's list + assert [t.name for t in binding.get_targets()] == ['width', 'area'] + + def test_diffusion_dict_targets_map_dataset_keys(self): + # WHEN: mapping DeltaLorentz predictions to custom dataset keys + model = DeltaLorentz() + binding = FitBinding( + model=model, + targets={ + 'width': 'L width', + 'area': 'L area', + 'delta_area': 'Elastic area', + }, + ) - def test_get_model_names_diffusion(self, diffusion_fit_binding): - # WHEN THEN - model_names = diffusion_fit_binding.get_model_names() + # THEN + targets = {t.name: t for t in binding.get_targets()} # EXPECT - assert model_names == [ - 'BrownianTranslationalDiffusion area', - 'BrownianTranslationalDiffusion width', - ] + assert targets['width'].dataset_key == 'L width' + assert targets['area'].dataset_key == 'L area' + assert targets['delta_area'].dataset_key == 'Elastic area' - def test_get_parameter_names(self, fit_binding): - # WHEN THEN - parameter_names = fit_binding.get_parameter_names() + def test_delta_lorentz_default_targets(self): + # WHEN: default keys derive from the component names + model = DeltaLorentz(lorentzian_name='Lorentzian', delta_name='Delta function') + binding = FitBinding(model=model) - # EXPECT - assert parameter_names == ['parameter1'] - - def test_get_parameter_names_diffusion(self, diffusion_fit_binding): - # WHEN THEN - parameter_names = diffusion_fit_binding.get_parameter_names() + # THEN + targets = {t.name: t for t in binding.get_targets()} # EXPECT - assert parameter_names == ['parameter3 area', 'parameter3 width'] + assert set(targets) == {'area', 'width', 'delta_area'} + assert targets['area'].dataset_key == 'Lorentzian area' + assert targets['width'].dataset_key == 'Lorentzian width' + assert targets['delta_area'].dataset_key == 'Delta function area' - # ------------------------------------------------------------------ - # Private methods - # ------------------------------------------------------------------ + def test_diffusion_target_units(self, diffusion_binding): + # WHEN + targets = {t.name: t for t in diffusion_binding.get_targets()} - def test_build_diffusion_callable(self, diffusion_fit_binding): + # THEN EXPECT: predictions take Q in 1/angstrom; widths are in the model's x_unit and + # areas in the scale unit + assert targets['width'].x_unit == '1/angstrom' + assert targets['width'].y_unit == diffusion_binding.model.x_unit + assert targets['area'].y_unit == str(diffusion_binding.model.scale.unit) - # WHEN - mock_model = Mock() - mock_model.calculate_QISF.return_value = 2.0 - mock_model.scale.value = 3.0 - mock_model.calculate_width.return_value = 0.5 - diffusion_fit_binding._model = mock_model + def test_diffusion_target_units_track_conversions(self, diffusion_binding): + # WHEN: targets are snapshots of the live model, so unit conversions are reflected + diffusion_binding.model.convert_x_unit('ueV') # THEN - area_callable = diffusion_fit_binding._build_diffusion_callable(mode='area') - width_callable = diffusion_fit_binding._build_diffusion_callable(mode='width') + targets = {t.name: t for t in diffusion_binding.get_targets()} # EXPECT - assert area_callable(0) == pytest.approx(6.0) # 2.0 * 3.0 - mock_model.calculate_QISF.assert_called_once_with(0) + assert targets['width'].y_unit == 'ueV' - assert width_callable(0) == pytest.approx(0.5) - mock_model.calculate_width.assert_called_once_with(0) + def test_diffusion_target_functions(self): + # WHEN + model = BrownianTranslationalDiffusion(diffusion_coefficient=2.4e-9, scale=0.5) + binding = FitBinding(model=model) + Q = np.array([0.5, 1.0]) # THEN - result_area = area_callable(0, unused_arg=123) - result_width = width_callable(0, unused_arg=123) + targets = {t.name: t for t in binding.get_targets()} # EXPECT - assert result_area == pytest.approx(6.0) # Should ignore unused_arg - assert result_width == pytest.approx(0.5) # Should ignore unused_arg - - def test_build_diffusion_callable_errors(self, diffusion_fit_binding): - with pytest.raises(ValueError, match='Unknown diffusion mode: invalid_mode'): - diffusion_fit_binding._build_diffusion_callable(mode='invalid_mode') + np.testing.assert_allclose(targets['width'].function(Q), model.calculate_width(Q)) + np.testing.assert_allclose( + targets['area'].function(Q), model.calculate_QISF(Q) * model.scale.value + ) - def test_get_modes(self, diffusion_fit_binding): + def test_delta_lorentz_delta_area_function(self): # WHEN - modes = diffusion_fit_binding._get_modes() - - # EXPECT - assert modes == ['area', 'width'] + model = DeltaLorentz(A_0=0.5, mean_u_squared=0.3, scale=2.0, lorentzian_width=0.1) + binding = FitBinding(model=model, targets=['delta_area']) + Q = np.array([0.5, 1.0]) # THEN - diffusion_fit_binding.modes = 'area' - modes = diffusion_fit_binding._get_modes() + target = binding.get_targets()[0] + # EXPECT - assert modes == ['area'] + np.testing.assert_allclose(target.function(Q), model.calculate_EISF(Q) * model.scale.value) # ------------------------------------------------------------------ # dunder methods # ------------------------------------------------------------------ - def test_repr(self, fit_binding): - # WHEN - repr_str = repr(fit_binding) - # THEN EXPECT + def test_repr(self, diffusion_binding): + # WHEN THEN + repr_str = repr(diffusion_binding) + + # EXPECT assert 'FitBinding' in repr_str - assert "parameter_name='parameter1'" in repr_str - assert "model='Gaussian'" in repr_str - assert 'modes=None' in repr_str + assert 'model=' in repr_str + assert 'targets=' in repr_str + + +class TestFitBindingWorkflows: + """End-to-end regression tests for the standard ParameterAnalysis workflows.""" + + def test_polynomial_targets_gaussian_area(self): + # WHEN: fitting a Polynomial to a 'Gaussian area' dataset key + polynomial = Polynomial( + coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='meV', display_name='Line' + ) + binding = FitBinding(model=polynomial, targets='Gaussian area') + + # THEN + target = binding.get_targets()[0] + + # EXPECT + assert target.dataset_key == 'Gaussian area' + assert target.label == 'Line' diff --git a/tests/unit/easydynamics/analysis/test_parameter_analysis.py b/tests/unit/easydynamics/analysis/test_parameter_analysis.py index 14d9c58c..031f813c 100644 --- a/tests/unit/easydynamics/analysis/test_parameter_analysis.py +++ b/tests/unit/easydynamics/analysis/test_parameter_analysis.py @@ -17,12 +17,25 @@ from easydynamics.sample_model.diffusion_model.brownian_translational_diffusion import ( BrownianTranslationalDiffusion, ) +from easydynamics.utils.fit_target import FitTarget + + +def make_target(dataset_key, function, label, x_unit=None, y_unit=None, name='value'): + """Build a FitTarget for mocking FitBinding.get_targets in tests.""" + return FitTarget( + name=name, + dataset_key=dataset_key, + function=function, + label=label, + x_unit=x_unit, + y_unit=y_unit, + ) class TestParameterAnalysis: @pytest.fixture def dataset(self): - Q = sc.array(dims=['Q'], values=[0.1, 0.2]) + Q = sc.array(dims=['Q'], values=[0.1, 0.2], unit='1/angstrom') return sc.Dataset( data={ 'parameter1': sc.DataArray( @@ -63,11 +76,14 @@ def mock_model_dataset(self): @pytest.fixture def parameter_analysis(self, dataset): - model = Polynomial(coefficients=[1.0, 0.5]) + model = Polynomial(coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='meV') diffusion_model = BrownianTranslationalDiffusion() - fit_binding1 = FitBinding(parameter_name='parameter1', model=model) - fit_binding2 = FitBinding(parameter_name='parameter3', model=diffusion_model) + fit_binding1 = FitBinding(model=model, targets='parameter1') + fit_binding2 = FitBinding( + model=diffusion_model, + targets={'area': 'parameter3 area', 'width': 'parameter3 width'}, + ) return ParameterAnalysis(parameters=dataset, bindings=[fit_binding1, fit_binding2]) @@ -75,8 +91,11 @@ def test_initialization(self, parameter_analysis): # WHEN THEN EXPECT assert isinstance(parameter_analysis, ParameterAnalysis) assert len(parameter_analysis.bindings) == 2 - assert parameter_analysis.bindings[0].parameter_name == 'parameter1' - assert parameter_analysis.bindings[1].parameter_name == 'parameter3' + assert parameter_analysis.bindings[0].targets == 'parameter1' + assert parameter_analysis.bindings[1].targets == { + 'area': 'parameter3 area', + 'width': 'parameter3 width', + } def test_parameter_property(self, parameter_analysis): # WHEN @@ -133,7 +152,7 @@ def test_bindings_property(self, parameter_analysis): # WHEN model = Polynomial(coefficients=[2.0, 1.0]) - new_binding = FitBinding(parameter_name='parameter2', model=model) + new_binding = FitBinding(model=model, targets='parameter2') parameter_analysis.bindings = new_binding # THEN EXPECT @@ -162,7 +181,7 @@ def test_fit_no_parameters_raises(self, parameter_analysis): def test_fit_wrong_parameter_name_raises(self, parameter_analysis): # WHEN model = Polynomial(coefficients=[2.0, 1.0]) - incorrect_binding = FitBinding(parameter_name='nonexistent_parameter', model=model) + incorrect_binding = FitBinding(model=model, targets='nonexistent_parameter') parameter_analysis.bindings = incorrect_binding # THEN EXPECT @@ -172,6 +191,87 @@ def test_fit_wrong_parameter_name_raises(self, parameter_analysis): ): parameter_analysis.fit() + def test_fit_incompatible_x_unit_raises(self, dataset): + # WHEN: Polynomial has x_unit='meV' but Q coordinate has unit '1/angstrom' + model = Polynomial(coefficients=[1.0, 0.5], x_unit='meV', y_unit='meV') + binding = FitBinding(model=model, targets='parameter1') + pa = ParameterAnalysis(parameters=dataset, bindings=binding) + + # THEN EXPECT + with pytest.raises(Exception, match=r'Q coordinate unit .* is incompatible'): + pa.fit() + + def test_fit_incompatible_y_unit_raises(self, dataset): + # WHEN: Polynomial has y_unit='1/angstrom' but parameter1 has unit 'meV' + model = Polynomial(coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='1/angstrom') + binding = FitBinding(model=model, targets='parameter1') + pa = ParameterAnalysis(parameters=dataset, bindings=binding) + + # THEN EXPECT + with pytest.raises(Exception, match="Parameter 'parameter1' unit 'meV' is incompatible"): + pa.fit() + + def test_fit_converts_x_unit(self): + # WHEN: Q is in '1/m' but the model declares x_unit='1/angstrom'. + # 1 1/m = 1e-10 1/angstrom, so values [1e10, 2e10] 1/m → [1.0, 2.0] 1/angstrom. + Q = sc.array(dims=['Q'], values=[1e10, 2e10], unit='1/m') + dataset = sc.Dataset( + data={ + 'param': sc.DataArray( + data=sc.array(dims=['Q'], values=[1.0, 2.0], variances=[0.1, 0.2], unit='meV'), + coords={'Q': Q}, + ) + } + ) + model = Polynomial(coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='meV') + pa = ParameterAnalysis( + parameters=dataset, bindings=FitBinding(model=model, targets='param') + ) + + with patch('easydynamics.analysis.parameter_analysis.MultiFitter') as MockMultiFitter: + MockMultiFitter.return_value.fit.return_value = MagicMock() + # THEN + pa.fit() + x_passed = MockMultiFitter.return_value.fit.call_args.kwargs['x'][0] + + # EXPECT: x values converted from 1/m to 1/angstrom + np.testing.assert_allclose(x_passed, [1.0, 2.0]) + + def test_fit_converts_y_unit(self): + # WHEN: parameter is in 'eV' but model declares y_unit='meV'. + # 1 eV = 1000 meV, so values [0.001, 0.002] eV → [1.0, 2.0] meV. + # Weights (= 1/std) scale by 1/y_factor: sqrt(1e-6)=1e-3 eV^-1 → 1.0 meV^-1. + Q = sc.array(dims=['Q'], values=[0.1, 0.2], unit='1/angstrom') + dataset = sc.Dataset( + data={ + 'param': sc.DataArray( + data=sc.array( + dims=['Q'], + values=[0.001, 0.002], + variances=[1e-6, 4e-6], + unit='eV', + ), + coords={'Q': Q}, + ) + } + ) + model = Polynomial(coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='meV') + pa = ParameterAnalysis( + parameters=dataset, bindings=FitBinding(model=model, targets='param') + ) + + with patch('easydynamics.analysis.parameter_analysis.MultiFitter') as MockMultiFitter: + MockMultiFitter.return_value.fit.return_value = MagicMock() + # THEN + pa.fit() + kwargs = MockMultiFitter.return_value.fit.call_args.kwargs + y_passed = kwargs['y'][0] + w_passed = kwargs['weights'][0] + + # EXPECT: y values converted from eV to meV; weights scaled inversely + np.testing.assert_allclose(y_passed, [1.0, 2.0]) + np.testing.assert_allclose(w_passed, [1.0, 0.5]) + def test_fit_success(self, parameter_analysis): # WHEN mock_result = MagicMock() @@ -469,9 +569,11 @@ def test_calculate_model_dataset_missing_parameter_raises(self, parameter_analys binding = parameter_analysis.bindings[0] with ( - patch.object(binding, 'get_parameter_names', return_value=['missing_name']), - patch.object(binding, 'get_model_names', return_value=['model']), - patch.object(binding, 'build_callables', return_value=[lambda x: x]), + patch.object( + binding, + 'get_targets', + return_value=[make_target('missing_name', lambda x: x, 'model')], + ), pytest.raises(ValueError, match='not found in parameters Dataset'), ): parameter_analysis.calculate_model_dataset([binding]) @@ -483,10 +585,10 @@ def test_calculate_model_dataset_single_binding(self, parameter_analysis): # Mock callable to return predictable output mock_callable = MagicMock(return_value=np.array([10.0, 20.0])) - with ( - patch.object(binding, 'build_callables', return_value=[mock_callable]), - patch.object(binding, 'get_model_names', return_value=['model1']), - patch.object(binding, 'get_parameter_names', return_value=['parameter1']), + with patch.object( + binding, + 'get_targets', + return_value=[make_target('parameter1', mock_callable, 'model1')], ): # THEN result = parameter_analysis.calculate_model_dataset([binding]) @@ -516,12 +618,16 @@ def test_calculate_model_dataset_multiple_bindings(self, parameter_analysis): mock_callable2 = MagicMock(return_value=np.array([30.0, 40.0])) with ( - patch.object(binding1, 'build_callables', return_value=[mock_callable1]), - patch.object(binding1, 'get_model_names', return_value=['model1']), - patch.object(binding1, 'get_parameter_names', return_value=['parameter1']), - patch.object(binding2, 'build_callables', return_value=[mock_callable2]), - patch.object(binding2, 'get_model_names', return_value=['model2']), - patch.object(binding2, 'get_parameter_names', return_value=['parameter2']), + patch.object( + binding1, + 'get_targets', + return_value=[make_target('parameter1', mock_callable1, 'model1')], + ), + patch.object( + binding2, + 'get_targets', + return_value=[make_target('parameter2', mock_callable2, 'model2')], + ), ): # THEN result = parameter_analysis.calculate_model_dataset([binding1, binding2]) @@ -562,19 +668,18 @@ def test_calculate_model_dataset_multiple_bindings_diffusion(self, parameter_ana mock_callable2w = MagicMock(return_value=np.array([50.0, 60.0])) with ( - patch.object(binding1, 'build_callables', return_value=[mock_callable1]), - patch.object(binding1, 'get_model_names', return_value=['model1']), - patch.object(binding1, 'get_parameter_names', return_value=['parameter1']), patch.object( - binding2, - 'build_callables', - return_value=[mock_callable2a, mock_callable2w], + binding1, + 'get_targets', + return_value=[make_target('parameter1', mock_callable1, 'model1')], ), - patch.object(binding2, 'get_model_names', return_value=['model2a', 'model2w']), patch.object( binding2, - 'get_parameter_names', - return_value=['parameter3 area', 'parameter3 width'], + 'get_targets', + return_value=[ + make_target('parameter3 area', mock_callable2a, 'model2a', name='area'), + make_target('parameter3 width', mock_callable2w, 'model2w', name='width'), + ], ), ): # THEN @@ -614,10 +719,76 @@ def test_calculate_model_dataset_multiple_bindings_diffusion(self, parameter_ana np.testing.assert_allclose(args[0], [0.1, 0.2]) assert kwargs == {} + def test_calculate_model_dataset_converts_x_unit(self): + # WHEN: Q is in '1/m' but model declares x_unit='1/angstrom'. + # Values [1e10, 2e10] 1/m → [1.0, 2.0] 1/angstrom passed to callable. + Q = sc.array(dims=['Q'], values=[1e10, 2e10], unit='1/m') + dataset = sc.Dataset( + data={ + 'param': sc.DataArray( + data=sc.array(dims=['Q'], values=[1.0, 2.0], unit='meV'), + coords={'Q': Q}, + ) + } + ) + model = Polynomial(coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='meV') + binding = FitBinding(model=model, targets='param') + pa = ParameterAnalysis(parameters=dataset, bindings=binding) + + mock_callable = MagicMock(return_value=np.array([10.0, 20.0])) + # THEN + with patch.object( + binding, + 'get_targets', + return_value=[ + make_target( + 'param', mock_callable, 'Polynomial', x_unit='1/angstrom', y_unit='meV' + ) + ], + ): + pa.calculate_model_dataset([binding]) + + # EXPECT: callable received x in 1/angstrom, not raw 1/m values + args, _ = mock_callable.call_args + np.testing.assert_allclose(args[0], [1.0, 2.0]) + + def test_calculate_model_dataset_converts_y_unit(self): + # WHEN: parameter is in 'eV' but model declares y_unit='meV'. + # Callable returns [1.0, 2.0] meV → stored as [0.001, 0.002] eV in the DataArray. + Q = sc.array(dims=['Q'], values=[0.1, 0.2], unit='1/angstrom') + dataset = sc.Dataset( + data={ + 'param': sc.DataArray( + data=sc.array(dims=['Q'], values=[0.001, 0.002], unit='eV'), + coords={'Q': Q}, + ) + } + ) + model = Polynomial(coefficients=[1.0, 0.5], x_unit='1/angstrom', y_unit='meV') + binding = FitBinding(model=model, targets='param') + pa = ParameterAnalysis(parameters=dataset, bindings=binding) + + mock_callable = MagicMock(return_value=np.array([1.0, 2.0])) # meV + # THEN + with patch.object( + binding, + 'get_targets', + return_value=[ + make_target( + 'param', mock_callable, 'Polynomial', x_unit='1/angstrom', y_unit='meV' + ) + ], + ): + result = pa.calculate_model_dataset([binding]) + + # EXPECT: model output converted from meV back to eV + np.testing.assert_allclose(result['Polynomial'].values, [0.001, 0.002]) + assert result['Polynomial'].unit == sc.Unit('eV') + def test_append_binding(self, parameter_analysis): # WHEN model = Polynomial(coefficients=[2.0, 1.0]) - new_binding = FitBinding(parameter_name='parameter2', model=model) + new_binding = FitBinding(model=model, targets='parameter2') # THEN parameter_analysis.append_binding(new_binding) @@ -666,8 +837,8 @@ def test_get_all_variables_overlapping_variables(self, parameter_analysis): # Create two bindings with overlapping variables model1 = Gaussian(display_name='model1') - binding1 = FitBinding(parameter_name='parameter1', model=model1) - binding2 = FitBinding(parameter_name='parameter2', model=model1) + binding1 = FitBinding(model=model1, targets='parameter1') + binding2 = FitBinding(model=model1, targets='parameter2') parameter_analysis.bindings = [binding1, binding2] @@ -832,14 +1003,14 @@ def test_get_xyweight_from_dataset_wrong_parameter_name_raises(self, parameter_a parameter_analysis._get_xyweight_from_dataset('nonexistent_parameter') @pytest.mark.parametrize( - 'non_finite_variance', - [np.inf, -np.inf, np.nan, -1.0, 0.0], - ids=['inf', '-inf', 'nan', 'negative', 'zero'], + 'bad_variance', + [np.inf, -np.inf, -1.0, 0.0], + ids=['inf', '-inf', 'negative', 'zero'], ) def test_get_xyweight_from_dataset_non_finite_weights_raises( - self, parameter_analysis, non_finite_variance + self, parameter_analysis, bad_variance ): - # WHEN + # Non-NaN invalid variances (inf, negative, zero) should still raise. Q = sc.array(dims=['Q'], values=[0.1, 0.2]) parameter_analysis.parameters = sc.Dataset( data={ @@ -847,7 +1018,7 @@ def test_get_xyweight_from_dataset_non_finite_weights_raises( data=sc.array( dims=['Q'], values=[1.0, 2.0], - variances=[1.0, non_finite_variance], + variances=[1.0, bad_variance], unit='meV', ), coords={'Q': Q}, @@ -861,6 +1032,33 @@ def test_get_xyweight_from_dataset_non_finite_weights_raises( ): parameter_analysis._get_xyweight_from_dataset('parameter1') + def test_get_xyweight_from_dataset_nan_variance_filters_row(self, parameter_analysis): + # NaN variances arise when a parameter is absent for a given Q; those rows are filtered. + + # WHEN + Q = sc.array(dims=['Q'], values=[0.1, 0.2]) + parameter_analysis.parameters = sc.Dataset( + data={ + 'parameter1': sc.DataArray( + data=sc.array( + dims=['Q'], + values=[1.0, np.nan], + variances=[0.25, np.nan], + unit='meV', + ), + coords={'Q': Q}, + ), + } + ) + + # THEN + x, y, w = parameter_analysis._get_xyweight_from_dataset('parameter1') + + # EXPECT + np.testing.assert_allclose(x, [0.1]) + np.testing.assert_allclose(y, [1.0]) + np.testing.assert_allclose(w, [1 / np.sqrt(0.25)]) + def test_get_xyweight_from_dataset_valid(self, parameter_analysis): # WHEN THEN x, y, w = parameter_analysis._get_xyweight_from_dataset('parameter1') @@ -871,6 +1069,74 @@ def test_get_xyweight_from_dataset_valid(self, parameter_analysis): expected_w = 1 / np.sqrt([0.1, 0.2]) np.testing.assert_allclose(w, expected_w) + def test_get_xyweight_from_dataset_no_variances(self, parameter_analysis): + # WHEN + Q = sc.array(dims=['Q'], values=[0.1, 0.2], unit='1/angstrom') + parameter_analysis.parameters = sc.Dataset( + data={ + 'parameter1': sc.DataArray( + data=sc.array(dims=['Q'], values=[1.0, 2.0], unit='meV'), + coords={'Q': Q}, + ), + } + ) + + # THEN + x, y, w = parameter_analysis._get_xyweight_from_dataset('parameter1') + + # EXPECT + np.testing.assert_allclose(x, [0.1, 0.2]) + np.testing.assert_allclose(y, [1.0, 2.0]) + np.testing.assert_allclose(w, [1.0, 1.0]) + + def test_get_xyweight_from_dataset_all_nan_variances_raises(self, parameter_analysis): + # WHEN + Q = sc.array(dims=['Q'], values=[0.1, 0.2], unit='1/angstrom') + parameter_analysis.parameters = sc.Dataset( + data={ + 'parameter1': sc.DataArray( + data=sc.array( + dims=['Q'], + values=[np.nan, np.nan], + variances=[np.nan, np.nan], + unit='meV', + ), + coords={'Q': Q}, + ), + } + ) + + # THEN EXPECT + with pytest.raises( + ValueError, + match="No finite positive variances found for parameter 'parameter1'", + ): + parameter_analysis._get_xyweight_from_dataset('parameter1') + + def test_get_unit_conversions_none_x_unit(self, parameter_analysis): + # WHEN + model = Polynomial(coefficients=[1.0], x_unit=None, y_unit='meV') + binding = FitBinding(model=model, targets='parameter1') + + # THEN + x_factor, y_factor = parameter_analysis._get_unit_conversions(binding.get_targets()[0]) + + # EXPECT + assert x_factor == pytest.approx(1.0) + assert y_factor == pytest.approx(1.0) + + def test_get_unit_conversions_none_y_unit(self, parameter_analysis): + # WHEN + model = Polynomial(coefficients=[1.0], x_unit='1/angstrom', y_unit=None) + binding = FitBinding(model=model, targets='parameter1') + + # THEN + x_factor, y_factor = parameter_analysis._get_unit_conversions(binding.get_targets()[0]) + + # EXPECT + assert x_factor == pytest.approx(1.0) + assert y_factor == pytest.approx(1.0) + def test_repr(self, parameter_analysis): # WHEN repr_str = repr(parameter_analysis) @@ -883,3 +1149,102 @@ def test_repr(self, parameter_analysis): assert f'n_parameters={len(parameter_analysis.parameters)}' in repr_str assert 'parameter_names=' in repr_str assert 'bindings=' in repr_str + + +class TestParameterAnalysisWorkflows: + """End-to-end fits for the standard workflows on synthetic data.""" + + @staticmethod + def _dataset_from_targets(model, Q, unit_overrides=None): + """Build a parameters Dataset from a model's own predictions (no noise).""" + unit_overrides = unit_overrides or {} + Q_coord = sc.array(dims=['Q'], values=Q, unit='1/angstrom') + data = {} + for target in model.get_fit_targets(): + values = target.function(Q) + unit = unit_overrides.get(target.name, target.y_unit) + factor = sc.to_unit(sc.scalar(1.0, unit=target.y_unit), unit).value + data[target.dataset_key] = sc.DataArray( + data=sc.array( + dims=['Q'], + values=values * factor, + variances=(0.01 * values * factor) ** 2, + unit=unit, + ), + coords={'Q': Q_coord}, + ) + return sc.Dataset(data) + + def test_delta_lorentz_three_target_simultaneous_fit(self): + # WHEN: synthetic width, area, and delta area curves from a known DeltaLorentz + from easydynamics.sample_model.diffusion_model.delta_lorentz import DeltaLorentz + + Q = np.linspace(0.4, 2.0, 9) + truth = DeltaLorentz(scale=2.0, mean_u_squared=0.3, A_0=0.6, lorentzian_width=0.12) + dataset = self._dataset_from_targets(truth, Q) + + # THEN: fit a fresh DeltaLorentz to all three dataset keys simultaneously + fit_model = DeltaLorentz(scale=1.5, mean_u_squared=0.2, A_0=0.5, lorentzian_width=0.1) + pa = ParameterAnalysis(parameters=dataset, bindings=FitBinding(model=fit_model)) + pa.fit() + + # EXPECT: the shared parameters are recovered across all three curves + assert fit_model.scale.value == pytest.approx(2.0, rel=1e-3) + assert fit_model.mean_u_squared.value == pytest.approx(0.3, rel=1e-3) + assert fit_model.A_0.value == pytest.approx(0.6, rel=1e-3) + assert fit_model.lorentzian_width.value == pytest.approx(0.12, rel=1e-3) + + def test_jump_diffusion_width_only_fit(self): + # WHEN: synthetic widths from a known jump diffusion model + from easydynamics.sample_model.diffusion_model.jump_translational_diffusion import ( + JumpTranslationalDiffusion, + ) + + Q = np.linspace(0.4, 2.0, 9) + truth = JumpTranslationalDiffusion(diffusion_coefficient=2.4e-9, relaxation_time=2.0) + dataset = self._dataset_from_targets(truth, Q) + + # THEN: fit only the width prediction + fit_model = JumpTranslationalDiffusion(diffusion_coefficient=1.0e-9, relaxation_time=1.0) + fit_model.scale.fixed = True + pa = ParameterAnalysis( + parameters=dataset, bindings=FitBinding(model=fit_model, targets=['width']) + ) + pa.fit() + + # EXPECT + assert fit_model.diffusion_coefficient.value == pytest.approx(2.4e-9, rel=1e-3) + assert fit_model.relaxation_time.value == pytest.approx(2.0, rel=1e-3) + + def test_diffusion_fit_converts_dataset_units(self): + # WHEN: the dataset stores widths in ueV and Q in 1/nm, while the model works in meV + # and 1/angstrom (regression: diffusion fits used to ignore units entirely, silently + # misfitting scaled data) + Q = np.linspace(0.4, 2.0, 9) + truth = BrownianTranslationalDiffusion(diffusion_coefficient=2.4e-9) + width_target = {t.name: t for t in truth.get_fit_targets()}['width'] + widths_mev = width_target.function(Q) + + Q_coord = sc.array(dims=['Q'], values=Q * 10, unit='1/nm') # 1 1/angstrom = 10 1/nm + dataset = sc.Dataset({ + width_target.dataset_key: sc.DataArray( + data=sc.array( + dims=['Q'], + values=widths_mev * 1000, # meV -> ueV + variances=(0.01 * widths_mev * 1000) ** 2, + unit='ueV', + ), + coords={'Q': Q_coord}, + ) + }) + + # THEN + fit_model = BrownianTranslationalDiffusion(diffusion_coefficient=1.0e-9) + fit_model.scale.fixed = True + pa = ParameterAnalysis( + parameters=dataset, bindings=FitBinding(model=fit_model, targets=['width']) + ) + pa.fit() + + # EXPECT: the diffusion coefficient is recovered despite the unit differences + assert fit_model.diffusion_coefficient.value == pytest.approx(2.4e-9, rel=1e-3)