From c4ca25cf5e5a47141f5d47809a22efe25640e774 Mon Sep 17 00:00:00 2001 From: Daniel Weindl Date: Mon, 18 Dec 2023 19:58:08 +0100 Subject: [PATCH 1/2] Add debugging helper get_model_for_preeq Makes life easier when debugging preequilibration issues. --- python/sdist/amici/debugging/__init__.py | 45 ++++++++++++++++++++++++ python/tests/test_preequilibration.py | 32 +++++++++++++++++ 2 files changed, 77 insertions(+) create mode 100644 python/sdist/amici/debugging/__init__.py diff --git a/python/sdist/amici/debugging/__init__.py b/python/sdist/amici/debugging/__init__.py new file mode 100644 index 0000000000..81663d17b9 --- /dev/null +++ b/python/sdist/amici/debugging/__init__.py @@ -0,0 +1,45 @@ +"""Functions for debugging AMICI simulation failures.""" +import amici +import numpy as np + + +def get_model_for_preeq(model: amici.Model, edata: amici.ExpData): + """Get a model set-up to simulate the preequilibration condition as + specified in `edata`. + + Useful for analyzing simulation errors during preequilibration. + Simulating the returned model will reproduce the behavior of + simulation-based preequilibration. + + Note that for models with events, the simulation results may differ. + During preequilibration, event-handling is disabled. However, when + simulating the returned model, event-handling will be enabled. + For events triggered at fixed timepoints, this can be avoided by setting + :meth:`t0 ` to a timepoints after the last trigger + timepoint. + + :param model: + The model for which *edata* was generated. + :param edata: + The experimental data object with a preequilibration condition + specified. + :return: + A copy of *model* with the same parameters, initial states, ... + as *amici_model* for the preequilibration condition. + Output timepoints are set to ``[inf]`` and will have to be adjusted. + """ + model = model.clone() + model.setTimepoints([np.inf]) + model.setFixedParameters(edata.fixedParametersPreequilibration) + if edata.pscale: + model.setParameterScale(edata.pscale) + if edata.parameters: + model.setParameters(edata.parameters) + if edata.plist: + model.setParameterList(edata.plist) + model.setInitialStates(edata.x0) + # has to be set *after* parameter list/scale! + model.setInitialStateSensitivities(edata.sx0) + model.setT0(edata.tstart_) + + return model diff --git a/python/tests/test_preequilibration.py b/python/tests/test_preequilibration.py index a42bc6354d..0c92a7190a 100644 --- a/python/tests/test_preequilibration.py +++ b/python/tests/test_preequilibration.py @@ -5,6 +5,7 @@ import amici import numpy as np import pytest +from amici.debugging import get_model_for_preeq from numpy.testing import assert_allclose from test_pysb import get_data @@ -633,3 +634,34 @@ def test_simulation_errors(preeq_fixture): assert rdata._swigptr.messages[2].identifier == "OTHER" assert rdata._swigptr.messages[3].severity == amici.LogSeverity_debug assert rdata._swigptr.messages[3].identifier == "BACKTRACE" + + +def test_get_model_for_preeq(preeq_fixture): + ( + model, + solver, + edata, + edata_preeq, + edata_presim, + edata_sim, + pscales, + plists, + ) = preeq_fixture + model.setSteadyStateSensitivityMode( + amici.SteadyStateSensitivityMode.integrationOnly + ) + model_preeq = get_model_for_preeq(model, edata) + rdata1 = amici.runAmiciSimulation(model_preeq, solver) + rdata2 = amici.runAmiciSimulation(model, solver, edata_preeq) + assert_allclose( + rdata1.x, + rdata2.x, + atol=solver.getAbsoluteTolerance(), + rtol=solver.getRelativeTolerance(), + ) + assert_allclose( + rdata1.sx, + rdata2.sx, + atol=solver.getAbsoluteTolerance(), + rtol=solver.getRelativeTolerance(), + ) From e022d1a888f83d3bd7a83fa653f70d653b21f703 Mon Sep 17 00:00:00 2001 From: Daniel Weindl Date: Tue, 19 Dec 2023 08:15:38 +0100 Subject: [PATCH 2/2] tol --- python/tests/test_preequilibration.py | 11 ++++------- 1 file changed, 4 insertions(+), 7 deletions(-) diff --git a/python/tests/test_preequilibration.py b/python/tests/test_preequilibration.py index 0c92a7190a..be447b0c54 100644 --- a/python/tests/test_preequilibration.py +++ b/python/tests/test_preequilibration.py @@ -6,7 +6,7 @@ import numpy as np import pytest from amici.debugging import get_model_for_preeq -from numpy.testing import assert_allclose +from numpy.testing import assert_allclose, assert_equal from test_pysb import get_data @@ -651,17 +651,14 @@ def test_get_model_for_preeq(preeq_fixture): amici.SteadyStateSensitivityMode.integrationOnly ) model_preeq = get_model_for_preeq(model, edata) + # the exactly same settings are used, so results should match exactly rdata1 = amici.runAmiciSimulation(model_preeq, solver) rdata2 = amici.runAmiciSimulation(model, solver, edata_preeq) - assert_allclose( + assert_equal( rdata1.x, rdata2.x, - atol=solver.getAbsoluteTolerance(), - rtol=solver.getRelativeTolerance(), ) - assert_allclose( + assert_equal( rdata1.sx, rdata2.sx, - atol=solver.getAbsoluteTolerance(), - rtol=solver.getRelativeTolerance(), )