From b28179970cd7ab0d9b7e276c3d25ff7dae69e7cd Mon Sep 17 00:00:00 2001 From: Joris Snellenburg Date: Fri, 27 Aug 2021 01:35:39 +0200 Subject: [PATCH] =?UTF-8?q?=F0=9F=A9=B9=20Fix=20and=20re-enable=20IRF=20Di?= =?UTF-8?q?spersion=20Test=20(#786)?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit * Fix and re-enable IRF Dispersion Test Improve the tests in test_spectral_irf Added NoDispersion test Re-enable IrfDispersion tests on the CI (de-Volkswagen the tests) Improve assignment of center_dispersion to result dataset * Rename center_dispersion to irf_center_location in result The name irf_center_location in the results describes more clearly what the variable holds. Deprecated irf model item center_dispersion for center_dispersion_coefficients Deprecated irf model item width_dispersion for width_dispersion_coefficients --- .github/workflows/CI_CD_actions.yml | 2 +- .../builtin/io/yml/test/test_model_parser.py | 4 +- .../builtin/io/yml/test/test_model_spec.yml | 35 ++-- glotaran/builtin/megacomplexes/decay/irf.py | 28 +-- .../decay/test/test_spectral_irf.py | 163 ++++++++++++++---- glotaran/builtin/megacomplexes/decay/util.py | 13 +- .../deprecation/modules/builtin_io_yml.py | 21 +++ .../modules/test/test_builtin_io_yml.py | 28 +++ 8 files changed, 223 insertions(+), 71 deletions(-) diff --git a/.github/workflows/CI_CD_actions.yml b/.github/workflows/CI_CD_actions.yml index f7b17499a..f8123d88c 100644 --- a/.github/workflows/CI_CD_actions.yml +++ b/.github/workflows/CI_CD_actions.yml @@ -109,7 +109,7 @@ jobs: run: pip freeze - name: Run tests run: | - pytest --cov=./ --cov-report term --cov-report xml --cov-config pyproject.toml -k 'not IrfDispersion' glotaran + pytest --cov=./ --cov-report term --cov-report xml --cov-config pyproject.toml glotaran - name: Codecov Upload uses: codecov/codecov-action@v2 diff --git a/glotaran/builtin/io/yml/test/test_model_parser.py b/glotaran/builtin/io/yml/test/test_model_parser.py index 8fa942b55..be4ac0ce3 100644 --- a/glotaran/builtin/io/yml/test/test_model_parser.py +++ b/glotaran/builtin/io/yml/test/test_model_parser.py @@ -126,9 +126,9 @@ def test_irf(model): assert irf.width == want want = [3] if i == 1 else [5, 6] if i == 2: - assert irf.center_dispersion == want + assert irf.center_dispersion_coefficients == want want = [7, 8] - assert irf.width_dispersion == want + assert irf.width_dispersion_coefficients == want want = [9] assert irf.scale == want assert irf.normalize == (i == 1) diff --git a/glotaran/builtin/io/yml/test/test_model_spec.yml b/glotaran/builtin/io/yml/test/test_model_spec.yml index cf5dfce21..2c485fa82 100644 --- a/glotaran/builtin/io/yml/test/test_model_spec.yml +++ b/glotaran/builtin/io/yml/test/test_model_spec.yml @@ -1,6 +1,5 @@ default-megacomplex: decay - dataset: dataset1: megacomplex: [cmplx1] @@ -23,36 +22,36 @@ irf: irf2: type: spectral-gaussian center: [1, 2] - width: [3,4] + width: [3, 4] scale: [9] normalize: false backsweep: true backsweep_period: 55 dispersion_center: 55 - center_dispersion: [5,6] - width_dispersion: [7,8] + center_dispersion_coefficients: [5, 6] + width_dispersion_coefficients: [7, 8] model_dispersion_with_wavenumber: true initial_concentration: inputD1: - compartments: [s1,s2,s3] - parameters: [1,2,3] + compartments: [s1, s2, s3] + parameters: [1, 2, 3] inputD2: - compartments: [s1,s2,s3] - parameters: [1,2,3] + compartments: [s1, s2, s3] + parameters: [1, 2, 3] # Convention matrix notation column = source, row = target compartment # (2,1) means from 1 to 2 k_matrix: km1: matrix: - (s1, s1): '1' - (s2, s1): '2' - (s1, s2): '3' - (s3, s1): '4' - (s1, s3): '5' - (s4, s1): '6' - (s1, s4): '7' + (s1, s1): "1" + (s2, s1): "2" + (s1, s2): "3" + (s3, s1): "4" + (s1, s3): "5" + (s4, s1): "6" + (s1, s4): "7" shape: shape1: @@ -63,9 +62,9 @@ shape: megacomplex: cmplx1: - k_matrix: [km1] # A megacomplex has one or more k-matrices + k_matrix: [km1] # A megacomplex has one or more k-matrices cmplx2: - k_matrix: [km2] + k_matrix: [km2] cmplx3: type: "spectral" shape: @@ -97,7 +96,7 @@ relations: - source: s1 target: s2 parameter: 8 - interval: [[1,100], [2,200]] + interval: [[1, 100], [2, 200]] weights: - datasets: [d1, d2] diff --git a/glotaran/builtin/megacomplexes/decay/irf.py b/glotaran/builtin/megacomplexes/decay/irf.py index 7ddbc17f6..14d7d2d78 100644 --- a/glotaran/builtin/megacomplexes/decay/irf.py +++ b/glotaran/builtin/megacomplexes/decay/irf.py @@ -47,10 +47,10 @@ class IrfMultiGaussian: one or more center of the irf as parameter indices width: one or more widths of the gaussian as parameter index - center_dispersion: + center_dispersion_coefficients: polynomial coefficients for the dispersion of the center as list of parameter indices. None for no dispersion. - width_dispersion: + width_dispersion_coefficients: polynomial coefficients for the dispersion of the width as parameter indices. None for no dispersion. @@ -124,8 +124,8 @@ class IrfGaussian(IrfMultiGaussian): @model_item( properties={ "dispersion_center": {"type": Parameter, "allow_none": True}, - "center_dispersion": {"type": List[Parameter], "default": []}, - "width_dispersion": {"type": List[Parameter], "default": []}, + "center_dispersion_coefficients": {"type": List[Parameter], "default": []}, + "width_dispersion_coefficients": {"type": List[Parameter], "default": []}, "model_dispersion_with_wavenumber": {"type": bool, "default": False}, }, has_type=True, @@ -149,12 +149,12 @@ class IrfSpectralMultiGaussian(IrfMultiGaussian): one or more center of the irf as parameter indices width: one or more widths of the gaussian as parameter index - center_dispersion: - polynomial coefficients for the dispersion of the - center as list of parameter indices. None for no dispersion. - width_dispersion: - polynomial coefficients for the dispersion of the - width as parameter indices. None for no dispersion. + center_dispersion_coefficients: + list of parameters with polynomial coefficients describing + the dispersion of the irf center location. None for no dispersion. + width_dispersion_coefficients: + list of parameters with polynomial coefficients describing + the dispersion of the width of the irf. None for no dispersion. """ @@ -173,16 +173,16 @@ def parameter(self, global_index: int, global_axis: np.ndarray): else (index - self.dispersion_center) / 100 ) - if len(self.center_dispersion) != 0: + if len(self.center_dispersion_coefficients) != 0: if self.dispersion_center is None: raise ModelError(f"No dispersion center defined for irf '{self.label}'") - for i, disp in enumerate(self.center_dispersion): + for i, disp in enumerate(self.center_dispersion_coefficients): centers += disp * np.power(dist, i + 1) - if len(self.width_dispersion) != 0: + if len(self.width_dispersion_coefficients) != 0: if self.dispersion_center is None: raise ModelError(f"No dispersion center defined for irf '{self.label}'") - for i, disp in enumerate(self.width_dispersion): + for i, disp in enumerate(self.width_dispersion_coefficients): widths = widths + disp * np.power(dist, i + 1) return centers, widths, scale, shift, backsweep, backsweep_period diff --git a/glotaran/builtin/megacomplexes/decay/test/test_spectral_irf.py b/glotaran/builtin/megacomplexes/decay/test/test_spectral_irf.py index 4256cbfa2..d75d636c2 100644 --- a/glotaran/builtin/megacomplexes/decay/test/test_spectral_irf.py +++ b/glotaran/builtin/megacomplexes/decay/test/test_spectral_irf.py @@ -1,4 +1,6 @@ +import warnings from copy import deepcopy +from textwrap import dedent import numpy as np import pytest @@ -35,6 +37,14 @@ sh1: type: one """ +MODEL_NO_IRF_DISPERSION = f"""\ +{MODEL_BASE} +irf: + irf1: + type: spectral-gaussian + center: irf.center + width: irf.width +""" MODEL_SIMPLE_IRF_DISPERSION = f"""\ {MODEL_BASE} irf: @@ -43,7 +53,7 @@ center: irf.center width: irf.width dispersion_center: irf.dispersion_center - center_dispersion: [irf.center_dispersion] + center_dispersion_coefficients: [irf.cdc1] """ MODEL_MULTI_IRF_DISPERSION = f"""\ {MODEL_BASE} @@ -53,8 +63,19 @@ center: [irf.center] width: [irf.width] dispersion_center: irf.dispersion_center - center_dispersion: [irf.center_dispersion1, irf.center_dispersion2] - width_dispersion: [irf.width_dispersion] + center_dispersion_coefficients: [irf.cdc1, irf.cdc2] + width_dispersion_coefficients: [irf.wdc1] +""" + +MODEL_MULTIPULSE_IRF_DISPERSION = f"""\ +{MODEL_BASE} +irf: + irf1: + type: spectral-multi-gaussian + center: [irf.center1, irf.center2] + width: [irf.width] + dispersion_center: irf.dispersion_center + center_dispersion_coefficients: [irf.cdc1, irf.cdc2, irf.cdc3] """ PARAMETERS_BASE = """\ @@ -64,62 +85,111 @@ - ['1', 0.5, {'non-negative': False}] """ +PARAMETERS_NO_IRF_DISPERSION = f"""\ +{PARAMETERS_BASE} +irf: + - ['center', 0.3] + - ['width', 0.1] +""" + PARAMETERS_SIMPLE_IRF_DISPERSION = f"""\ {PARAMETERS_BASE} irf: - ['center', 0.3] - ['width', 0.1] - ['dispersion_center', 400, {{'vary': False}}] - - ['center_dispersion', 0.5] + - ['cdc1', 0.5] """ +# What is this? PARAMETERS_MULTI_IRF_DISPERSION = f"""\ {PARAMETERS_BASE} irf: - ["center", 0.3] - ["width", 0.1] - ["dispersion_center", 400, {{"vary": False}}] - - ["center_dispersion1", 0.01] - - ["center_dispersion2", 0.001] - - ["width_dispersion", 0.025] + - ["cdc1", 0.1] + - ["cdc2", 0.01] + - ["wdc1", 0.025] """ +PARAMETERS_MULTIPULSE_IRF_DISPERSION = f"""\ +{PARAMETERS_BASE} +irf: + - ["center1", 0.3] + - ["center2", 0.4] + - ['width', 0.1] + - ['dispersion_center', 400, {{'vary': False}}] + - ["cdc1", 0.5] + - ["cdc2", 0.1] + - ["cdc3", -0.01] +""" + + +def _time_axis(): + time_p1 = np.linspace(-1, 1, 20, endpoint=False) + time_p2 = np.linspace(1, 2, 10, endpoint=False) + time_p3 = np.geomspace(2, 20, num=20) + return np.array(np.concatenate([time_p1, time_p2, time_p3])) + + +def _spectral_axis(): + return np.linspace(300, 500, 3) + + +def _calculate_irf_position( + index, center, dispersion_center=None, center_dispersion_coefficients=None +): + if center_dispersion_coefficients is None: + center_dispersion_coefficients = [] + if dispersion_center is not None: + distance = (index - dispersion_center) / 100 + if dispersion_center is not None: + for i, coefficient in enumerate(center_dispersion_coefficients): + center += coefficient * np.power(distance, i + 1) + return center + + +class NoIrfDispersion: + model = load_model(MODEL_NO_IRF_DISPERSION, format_name="yml_str") + parameters = load_parameters(PARAMETERS_NO_IRF_DISPERSION, format_name="yml_str") + axis = {"time": _time_axis(), "spectral": _spectral_axis()} + class SimpleIrfDispersion: model = load_model(MODEL_SIMPLE_IRF_DISPERSION, format_name="yml_str") parameters = load_parameters(PARAMETERS_SIMPLE_IRF_DISPERSION, format_name="yml_str") - time_p1 = np.linspace(-1, 2, 50, endpoint=False) - time_p2 = np.linspace(2, 5, 30, endpoint=False) - time_p3 = np.geomspace(5, 10, num=20) - time = np.array(np.concatenate([time_p1, time_p2, time_p3])) - spectral = np.arange(300, 500, 100) - axis = {"time": time, "spectral": spectral} + axis = {"time": _time_axis(), "spectral": _spectral_axis()} class MultiIrfDispersion: model = load_model(MODEL_MULTI_IRF_DISPERSION, format_name="yml_str") parameters = load_parameters(PARAMETERS_MULTI_IRF_DISPERSION, format_name="yml_str") - time = np.arange(-1, 5, 0.2) - spectral = np.arange(300, 500, 100) - axis = {"time": time, "spectral": spectral} + axis = {"time": _time_axis(), "spectral": _spectral_axis()} + + +class MultiCenterIrfDispersion: + model = load_model(MODEL_MULTIPULSE_IRF_DISPERSION, format_name="yml_str") + parameters = load_parameters(PARAMETERS_MULTIPULSE_IRF_DISPERSION, format_name="yml_str") + axis = {"time": _time_axis(), "spectral": _spectral_axis()} @pytest.mark.parametrize( "suite", [ + NoIrfDispersion, SimpleIrfDispersion, MultiIrfDispersion, + MultiCenterIrfDispersion, ], ) def test_spectral_irf(suite): model = suite.model - print(model.validate()) - assert model.valid() + assert model.valid(), model.validate() parameters = suite.parameters - print(model.validate(parameters)) - assert model.valid(parameters) + assert model.valid(parameters), model.validate(parameters) sim_model = deepcopy(model) sim_model.dataset["dataset1"].global_megacomplex = ["mc2"] @@ -136,27 +206,62 @@ def test_spectral_irf(suite): maximum_number_function_evaluations=20, ) result = optimize(scheme) - print(result.optimized_parameters) + # print(result.optimized_parameters) for label, param in result.optimized_parameters.all(): assert np.allclose(param.value, parameters.get(label).value, rtol=1e-1) resultdata = result.data["dataset1"] - print(resultdata) - + # print(resultdata) assert np.array_equal(dataset["time"], resultdata["time"]) assert np.array_equal(dataset["spectral"], resultdata["spectral"]) assert dataset.data.shape == resultdata.data.shape assert dataset.data.shape == resultdata.fitted_data.shape - assert np.allclose(dataset.data, resultdata.fitted_data, atol=1e-14) + # assert np.allclose(dataset.data, resultdata.fitted_data, atol=1e-14) + + fit_data_max_at_start = resultdata.fitted_data.isel(spectral=0).argmax(axis=0) + fit_data_max_at_end = resultdata.fitted_data.isel(spectral=-1).argmax(axis=0) + + if suite is NoIrfDispersion: + assert "center_dispersion_1" not in resultdata + assert fit_data_max_at_start == fit_data_max_at_end + else: + assert "center_dispersion_1" in resultdata + assert fit_data_max_at_start != fit_data_max_at_end + if abs(fit_data_max_at_start - fit_data_max_at_end) < 3: + warnings.warn( + dedent( + """ + Bad test, one of the following could be the case: + - dispersion too small + - spectral window to small + - time resolution (around the maximum of the IRF) too low" + """ + ) + ) - irf_max_at_start = resultdata.fitted_data.isel(spectral=0).argmax(axis=0) - irf_max_at_end = resultdata.fitted_data.isel(spectral=-1).argmax(axis=0) - print(f" irf_max_at_start: {irf_max_at_start}\n irf_max_at_end: {irf_max_at_end}") - # These should not be equal due to dispersion: - assert irf_max_at_start != irf_max_at_end + for x in suite.axis["spectral"]: + # calculated irf location + model_irf_center = suite.model.irf["irf1"].center + model_dispersion_center = suite.model.irf["irf1"].dispersion_center + model_center_dispersion_coefficients = suite.model.irf[ + "irf1" + ].center_dispersion_coefficients + calc_irf_location_at_x = _calculate_irf_position( + x, model_irf_center, model_dispersion_center, model_center_dispersion_coefficients + ) + # fitted irf location + fitted_irf_loc_at_x = resultdata["irf_center_location"].sel(spectral=x) + assert np.allclose(calc_irf_location_at_x, fitted_irf_loc_at_x) assert "species_associated_spectra" in resultdata assert "decay_associated_spectra" in resultdata assert "irf_center" in resultdata + + +if __name__ == "__main__": + test_spectral_irf(NoIrfDispersion) + test_spectral_irf(SimpleIrfDispersion) + test_spectral_irf(MultiIrfDispersion) + test_spectral_irf(MultiCenterIrfDispersion) diff --git a/glotaran/builtin/megacomplexes/decay/util.py b/glotaran/builtin/megacomplexes/decay/util.py index c51a3ce67..869d9d39c 100644 --- a/glotaran/builtin/megacomplexes/decay/util.py +++ b/glotaran/builtin/megacomplexes/decay/util.py @@ -228,10 +228,9 @@ def retrieve_irf(dataset_model: DatasetModel, dataset: xr.Dataset, global_dimens dataset["irf_center"] = ("irf_nr", center) if len(center) > 1 else center[0] dataset["irf_width"] = ("irf_nr", width) if len(width) > 1 else width[0] if isinstance(irf, IrfSpectralMultiGaussian) and irf.dispersion_center: - for i, dispersion in enumerate( - irf.calculate_dispersion(dataset.coords["spectral"].values) - ): - dataset[f"center_dispersion_{i+1}"] = ( - global_dimension, - dispersion, - ) + dataset["irf_center_location"] = ( + ("irf_nr", global_dimension), + irf.calculate_dispersion(dataset.coords["spectral"].values), + ) + # center_dispersion_1 for backwards compatibility (0.3-0.4.1) + dataset["center_dispersion_1"] = dataset["irf_center_location"].sel(irf_nr=0) diff --git a/glotaran/deprecation/modules/builtin_io_yml.py b/glotaran/deprecation/modules/builtin_io_yml.py index 635a5a9a9..0a3aa99dc 100644 --- a/glotaran/deprecation/modules/builtin_io_yml.py +++ b/glotaran/deprecation/modules/builtin_io_yml.py @@ -85,3 +85,24 @@ def model_spec_deprecations(spec: MutableMapping[Any, Any]) -> None: swap_keys=("equal_area_penalties", "clp_area_penalties"), stacklevel=load_model_stack_level, ) + + if "irf" in spec: + for _, irf in spec["irf"].items(): + deprecate_dict_entry( + dict_to_check=irf, + deprecated_usage="center_dispersion", + new_usage="center_dispersion_coefficients", + to_be_removed_in_version="0.7.0", + swap_keys=("center_dispersion", "center_dispersion_coefficients"), + stacklevel=load_model_stack_level, + ) + + for _, irf in spec["irf"].items(): + deprecate_dict_entry( + dict_to_check=irf, + deprecated_usage="width_dispersion", + new_usage="width_dispersion_coefficients", + to_be_removed_in_version="0.7.0", + swap_keys=("width_dispersion", "width_dispersion_coefficients"), + stacklevel=load_model_stack_level, + ) diff --git a/glotaran/deprecation/modules/test/test_builtin_io_yml.py b/glotaran/deprecation/modules/test/test_builtin_io_yml.py index 26202b947..8cf03f510 100644 --- a/glotaran/deprecation/modules/test/test_builtin_io_yml.py +++ b/glotaran/deprecation/modules/test/test_builtin_io_yml.py @@ -55,6 +55,32 @@ "clp_area_penalties", [{"type": "equal_area"}], ), + ( + dedent( + """ + irf: + irf1: + center_dispersion: [cdc1] + + """ + ), + 1, + "irf", + {"irf1": {"center_dispersion_coefficients": ["cdc1"]}}, + ), + ( + dedent( + """ + irf: + irf1: + "width_dispersion": [wdc1] + + """ + ), + 1, + "irf", + {"irf1": {"width_dispersion_coefficients": ["wdc1"]}}, + ), ), ids=( "type: kinetic-spectrum", @@ -62,6 +88,8 @@ "spectral_relations", "spectral_constraints", "equal_area_penalties", + "center_dispersion", + "width_dispersion", ), ) def test_model_spec_deprecations(