From c0d7a4caba6de4a2cbf25effad2d9623775def2e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Fri, 16 Oct 2026 09:10:16 +0100 Subject: [PATCH 01/15] assignments_to_experiments: convert values that ramp with time, via a rate rule A condition change is applied once at the start of its period, with `time` set to the period's start, so a value that changes with time within a period cannot be a condition change on a constant target, and such assignments were rejected. Fujita_SciSignal2010's EGF input is one: `piecewise(EGF_rate * time + EGF_0, time <= EGF_end, 0)`, steps in training and ramps in its held-out data. A value may now change linearly with time. The target then becomes a differential entity, with the rate rule `d target/dt = _slope` and a new constant parameter for the slope, and each period sets both the target's value at the period start and its slope. During preequilibration the value is held at its value at time 0 with slope zero, and the first period after preequilibration sets the slope again. An assignment that only switches between constant values converts exactly as before, as a constant parameter. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01MNJePorYs9zmPYMam6g4EP --- petab/v2/converters.py | 81 +++++++++++++++++--- tests/v2/test_conversion.py | 145 +++++++++++++++++++++++++++++++++++- 2 files changed, 213 insertions(+), 13 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 3c28a139..48ee2cdb 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -755,7 +755,8 @@ def assignments_to_experiments(problem: Problem) -> Problem: :data:`DEFAULT_EXPERIMENT_ID`, if needed. **Scope:** conversion is supported for SBML assignment rules that depend - on time only through the conditions of piecewise functions, where + on time through the conditions of piecewise functions, and at most + linearly through their values, where * every such condition is an inequality that can be re-arranged explicitly into `t />= ...`, @@ -770,7 +771,13 @@ def assignments_to_experiments(problem: Problem) -> Problem: given at different times (`d1 * piecewise(0, t < t1, 1) + d2 * piecewise(0, t < t2, 1)`), phases one after another as several pieces or nested piecewise functions (`piecewise(b0, t < t1, piecewise(b1, t < t2, - b2))`). Each switch that changes the value becomes a period. Anything else + b2))`). Each switch that changes the value becomes a period. A value may + also change linearly with time within a period (a ramp: `piecewise(r * t + + a, t <= t_end, 0)`). A condition change is applied once at the start of + its period, so the target of such an assignment becomes a differential + entity instead, with the rate rule `d target/dt = _slope`, and + each period sets both the target and its slope; an assignment that only + switches between constant values stays a constant parameter. Anything else raises a :class:`NotImplementedError`. .. warning:: @@ -855,7 +862,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: ) model_values = { s: sbml_model.getParameter(str(s)).getValue() - for s in default.free_symbols + for s in default.free_symbols - {time} } # unless a condition sets them, the switching times' symbols keep # their value from the parameter table or else from the model @@ -867,7 +874,24 @@ def assignments_to_experiments(problem: Problem) -> Problem: for s in t_switch.free_symbols } target = sbml_model.getParameter(target_id) - target.setValue(float(default.subs(model_values))) + # A value that changes with time within a period (a ramp) cannot be a + # condition change on a constant target, which is applied once at + # the start of its period. The target becomes a differential entity + # instead, `d target/dt = _slope`, and each period sets both. + ramps = expr.xreplace({c: sp.Dummy() for c in switches}).has(time) + slope_id = None + if ramps: + slope_id = f"{target_id}_slope" + while sbml_model.getElementBySId(slope_id) is not None: + slope_id += "_" + slope = sbml_model.createParameter() + slope.setId(slope_id) + slope.setConstant(False) + slope.setValue(float(sp.diff(default, time).subs(model_values))) + rate_rule = sbml_model.createRateRule() + rate_rule.setVariable(target_id) + rate_rule.setMath(libsbml.parseL3Formula(slope_id)) + target.setValue(float(default.subs(time, 0).subs(model_values))) target.setConstant(False) assignment.removeFromParentAndDelete() symbol_ids = {str(s) for s in expr.free_symbols} @@ -899,11 +923,16 @@ def assignments_to_experiments(problem: Problem) -> Problem: if period.is_preequilibration: # the value at time 0, which respects whether a - # condition holds for a switch at time 0 + # condition holds for a switch at time 0, held there: + # a ramp does not run during preequilibration value = expr.subs(subs).subs(time, 0) additions.append( (experiment, period.time, target_id, value) ) + if ramps: + additions.append( + (experiment, period.time, slope_id, sp.S.Zero) + ) continue t_end = periods[i + 1].time if i + 1 < len(periods) else sp.oo @@ -933,14 +962,27 @@ def assignments_to_experiments(problem: Problem) -> Problem: for _, value in switch_values ) ) + after_preequilibration = ( + i > 0 and periods[i - 1].is_preequilibration + ) if not redundant and ( - i == 0 or {c.target_id for c in changes} & symbol_ids + i == 0 + or {c.target_id for c in changes} & symbol_ids + # preequilibration held the slope at zero + or (ramps and after_preequilibration) ): - additions.append( - (experiment, period.time, target_id, start_value) + additions += _changes( + experiment, + period.time, + target_id, + start_value, + time, + slope_id, ) for t, value in switch_values: - additions.append((experiment, t, target_id, value)) + additions += _changes( + experiment, t, target_id, value, time, slope_id + ) if converted_ids and any(e.has_preequilibration for e in experiments): warnings.warn( @@ -1012,6 +1054,16 @@ def assignments_to_experiments(problem: Problem) -> Problem: return problem +def _changes(experiment, t, target_id, value, time, slope_id): + """The changes that start ``value`` at time ``t``: the target set to the + value at ``t``, and, for a target that ramps, its slope set to the + value's rate of change.""" + changes = [(experiment, t, target_id, value.subs(time, t))] + if slope_id is not None: + changes.append((experiment, t, slope_id, sp.diff(value, time))) + return changes + + def _value_from( expr: sp.Expr, switches: dict, @@ -1072,9 +1124,14 @@ def _supported_time_switches( raise unsupported t_switch = sp.solve(condition.gts - condition.lts, time)[0] switches[condition] = (t_switch, slope) - # time enters only through the conditions - if expr.xreplace({c: sp.Dummy() for c in switches}).has(time): - raise unsupported + # outside the conditions, time enters at most linearly: the value may + # ramp, at a rate that does not itself change with time + rest = expr.xreplace({c: sp.Dummy() for c in switches}) + if rest.has(time) and sp.diff(rest, time).has(time): + raise NotImplementedError( + "Only values that change at most linearly with time are " + f"supported, but the assignment to `{target_id}` is `{expr}`." + ) # no period can start at an estimated time if estimated := sorted( { diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index 7fd8fd96..ec9e7983 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -1,6 +1,7 @@ import logging from pathlib import Path +import libsbml import numpy as np import pandas as pd import pytest @@ -957,7 +958,8 @@ def test_assignments_to_experiments_preequilibration(): ("piecewise(0, time == 10, 1)", "single"), ("piecewise(0, time < 10 && time > 5, 1)", "single"), ("piecewise(0, time^2 < 10, 1)", "single"), - ("piecewise(time, time < 10, 1)", "single"), + ("piecewise(time^2, time < 10, 1)", "linearly"), + ("piecewise(exp(time), time < 10, 1)", "linearly"), ( "piecewise(0, time < t_dose, 1)", r"estimated parameters \['t_dose'\]", @@ -1196,6 +1198,147 @@ def test_assignments_to_experiments_condition_specific_doses(): } +def _rate_rule(problem, target_id): + rule = problem.model.sbml_model.getRateRuleByVariable(target_id) + return None if rule is None else libsbml.formulaToL3String(rule.getMath()) + + +def _changes_by_period(problem): + """{(experiment ID, time): {target ID: target value}}""" + conditions = problem.condition_df.set_index(v2.C.CONDITION_ID) + out = {} + for _, row in problem.experiment_df.iterrows(): + change = conditions.loc[[row[v2.C.CONDITION_ID]]] + out.setdefault((row[v2.C.EXPERIMENT_ID], row[v2.C.TIME]), {}).update( + { + r[v2.C.TARGET_ID]: float(r[v2.C.TARGET_VALUE]) + for _, r in change.iterrows() + } + ) + return out + + +def test_assignments_to_experiments_ramp(): + """Test that a value changing linearly with time becomes a rate rule + whose slope each period sets.""" + problem = Problem() + ant_model = """ + model conversion + species A, B; + A = 10; + B = 0; + r = 0.5; a = 1; + k1 := piecewise(r * time + a, time <= 10, 0); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement("obs_a", time=20.0, measurement=5.0) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + converted = assignments_to_experiments(problem) + + sbml = converted.model.sbml_model + # the ramp's start and slope are the model's own initial values + assert _rate_rule(converted, "k1") == "k1_slope" + assert sbml.getParameter("k1").getValue() == 1.0 + assert sbml.getParameter("k1_slope").getValue() == 0.5 + # at the switch, the target and its slope drop to zero + assert _changes_by_period(converted) == { + (DEFAULT_EXPERIMENT_ID, 10.0): {"k1": 0.0, "k1_slope": 0.0} + } + + +def test_assignments_to_experiments_condition_specific_ramps(): + """Test pulses and ramps whose start, rate and end are set by the + condition table, as in Fujita_SciSignal2010.""" + problem = Problem() + ant_model = """ + model conversion + species A, B; + A = 10; + B = 0; + EGF_0 = 0; EGF_rate = 0; EGF_end = 0; + k1 := piecewise(EGF_rate * time + EGF_0, time <= EGF_end, 0); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + # a pulse: a constant dose withdrawn at t = 60 + problem.add_condition("pulse", EGF_0=1.0, EGF_rate=0.0, EGF_end=60.0) + # a ramp: growing from zero at 0.01 per unit time until t = 3600 + problem.add_condition("ramp", EGF_0=0.0, EGF_rate=0.01, EGF_end=3600.0) + for c in ("pulse", "ramp"): + problem.add_experiment(f"experiment____{c}", 0, c) + problem.add_measurement( + "obs_a", + experiment_id=f"experiment____{c}", + time=30.0, + measurement=5.0, + ) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + converted = assignments_to_experiments(problem) + + assert _rate_rule(converted, "k1") == "k1_slope" + assert _changes_by_period(converted) == { + ("experiment____pulse", 0.0): {"k1": 1.0, "k1_slope": 0.0}, + ("experiment____pulse", 60.0): {"k1": 0.0, "k1_slope": 0.0}, + ("experiment____ramp", 0.0): {"k1": 0.0, "k1_slope": 0.01}, + ("experiment____ramp", 3600.0): {"k1": 0.0, "k1_slope": 0.0}, + } + + +def test_assignments_to_experiments_ramp_preequilibration(): + """Test that a ramp is held at its value at time 0 during + preequilibration and starts at time 0.""" + problem = Problem() + ant_model = """ + model conversion + species A, B; + A = 10; + B = 0; + r = 0; + k1 := piecewise(r * time + 1, time <= 10, 0); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + problem.add_condition("pre", r=0.0) + problem.add_condition("main", r=2.0) + problem.add_experiment( + "experiment__pre___main", v2.C.TIME_PREEQUILIBRATION, "pre", 0, "main" + ) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", + experiment_id="experiment__pre___main", + time=20.0, + measurement=5.0, + ) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + with pytest.warns(UserWarning, match="during preequilibration"): + converted = assignments_to_experiments(problem) + + assert _changes_by_period(converted) == { + ("experiment__pre___main", v2.C.TIME_PREEQUILIBRATION): { + "k1": 1.0, + "k1_slope": 0.0, + }, + ("experiment__pre___main", 0.0): {"k1": 1.0, "k1_slope": 2.0}, + ("experiment__pre___main", 10.0): {"k1": 0.0, "k1_slope": 0.0}, + } + + @pytest.mark.parametrize("t_switch", [0, 10]) @pytest.mark.parametrize( "condition, holds_at_switch, holds_before", From 368251b9270197c96cc16a4a26971cb67b747eb5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Fri, 16 Oct 2026 09:49:12 +0100 Subject: [PATCH 02/15] assignments_to_experiments: model change only makes a ramp's target differential The slope parameter is a placeholder at 0; every ramp's start value and slope, including at time 0, are set through the condition table. Constant model parameters in changes at an experiment's start take their model values, since those changes may only refer to parameter table symbols. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01MNJePorYs9zmPYMam6g4EP --- petab/v2/converters.py | 24 +++++++++++++++++++----- tests/v2/test_conversion.py | 12 +++++++----- 2 files changed, 26 insertions(+), 10 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 48ee2cdb..0ef9a14f 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -876,8 +876,9 @@ def assignments_to_experiments(problem: Problem) -> Problem: target = sbml_model.getParameter(target_id) # A value that changes with time within a period (a ramp) cannot be a # condition change on a constant target, which is applied once at - # the start of its period. The target becomes a differential entity - # instead, `d target/dt = _slope`, and each period sets both. + # the start of its period. The model only makes the target a + # differential entity, `d target/dt = _slope`; the slope's + # value, like the target's, is set by every period that starts one. ramps = expr.xreplace({c: sp.Dummy() for c in switches}).has(time) slope_id = None if ramps: @@ -887,7 +888,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: slope = sbml_model.createParameter() slope.setId(slope_id) slope.setConstant(False) - slope.setValue(float(sp.diff(default, time).subs(model_values))) + slope.setValue(0.0) rate_rule = sbml_model.createRateRule() rate_rule.setVariable(target_id) rate_rule.setMath(libsbml.parseL3Formula(slope_id)) @@ -951,9 +952,11 @@ def assignments_to_experiments(problem: Problem) -> Problem: # the first period, unless that value refers to parameter # table symbols, for which the model only has a number. The # switches can become the first periods if they only refer - # to parameter table symbols. + # to parameter table symbols. A ramp's start and slope are + # always set by the tables. redundant = ( - i == 0 + not ramps + and i == 0 and (start_value - default).is_zero and not {str(s) for s in default.free_symbols} & parameter_ids @@ -971,6 +974,17 @@ def assignments_to_experiments(problem: Problem) -> Problem: # preequilibration held the slope at zero or (ramps and after_preequilibration) ): + if i == 0: + # at the start of an experiment, a change may only + # refer to parameter table symbols; other symbols + # are constant model parameters, at their values + start_value = start_value.subs( + { + s: sbml_model.getParameter(str(s)).getValue() + for s in start_value.free_symbols - {time} + if str(s) not in parameter_ids + } + ) additions += _changes( experiment, period.time, diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index ec9e7983..306e63f9 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -1242,13 +1242,15 @@ def test_assignments_to_experiments_ramp(): converted = assignments_to_experiments(problem) sbml = converted.model.sbml_model - # the ramp's start and slope are the model's own initial values + # the model only makes the target differential; the slope is a + # placeholder until the tables set it assert _rate_rule(converted, "k1") == "k1_slope" - assert sbml.getParameter("k1").getValue() == 1.0 - assert sbml.getParameter("k1_slope").getValue() == 0.5 - # at the switch, the target and its slope drop to zero + assert sbml.getParameter("k1_slope").getValue() == 0.0 + # the ramp starts at time 0, and at the switch, the target and its slope + # drop to zero assert _changes_by_period(converted) == { - (DEFAULT_EXPERIMENT_ID, 10.0): {"k1": 0.0, "k1_slope": 0.0} + (DEFAULT_EXPERIMENT_ID, 0.0): {"k1": 1.0, "k1_slope": 0.5}, + (DEFAULT_EXPERIMENT_ID, 10.0): {"k1": 0.0, "k1_slope": 0.0}, } From 304dfe7a8764b57be236ccb3ff765b8ce8a9284b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Mon, 19 Oct 2026 19:32:12 +0100 Subject: [PATCH 03/15] assignments_to_experiments: rename the period index to `i_period` Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 0ef9a14f..7dcd5a66 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -902,7 +902,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: for experiment in experiments: periods = experiment.sorted_periods overrides = {} - for i, period in enumerate(periods): + for i_period, period in enumerate(periods): changes = problem.get_changes_for_period(period) overrides |= {c.target_id: c.target_value for c in changes} subs = fixed_values | { @@ -936,7 +936,11 @@ def assignments_to_experiments(problem: Problem) -> Problem: ) continue - t_end = periods[i + 1].time if i + 1 < len(periods) else sp.oo + t_end = ( + periods[i_period + 1].time + if i_period + 1 < len(periods) + else sp.oo + ) start_value = _value_from( expr, switches, time, subs, times, period.time ) @@ -956,7 +960,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: # always set by the tables. redundant = ( not ramps - and i == 0 + and i_period == 0 and (start_value - default).is_zero and not {str(s) for s in default.free_symbols} & parameter_ids @@ -966,15 +970,15 @@ def assignments_to_experiments(problem: Problem) -> Problem: ) ) after_preequilibration = ( - i > 0 and periods[i - 1].is_preequilibration + i_period > 0 and periods[i_period - 1].is_preequilibration ) if not redundant and ( - i == 0 + i_period == 0 or {c.target_id for c in changes} & symbol_ids # preequilibration held the slope at zero or (ramps and after_preequilibration) ): - if i == 0: + if i_period == 0: # at the start of an experiment, a change may only # refer to parameter table symbols; other symbols # are constant model parameters, at their values From 54328e6ad52254ad56c455d4fc063adc22d87d6c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 12:53:12 +0000 Subject: [PATCH 04/15] assignments_to_experiments: add a ramp's slope with `add_sbml_parameter` `add_sbml_parameter` checks libsbml's return codes, which the hand-written parameter creation did not. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 7dcd5a66..55bbdb92 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -885,10 +885,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: slope_id = f"{target_id}_slope" while sbml_model.getElementBySId(slope_id) is not None: slope_id += "_" - slope = sbml_model.createParameter() - slope.setId(slope_id) - slope.setConstant(False) - slope.setValue(0.0) + add_sbml_parameter(sbml_model, slope_id, value=0.0, constant=False) rate_rule = sbml_model.createRateRule() rate_rule.setVariable(target_id) rate_rule.setMath(libsbml.parseL3Formula(slope_id)) From b087dfc5f458e5929220852f26a395df412f3381 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 18 Oct 2026 22:38:12 +0100 Subject: [PATCH 05/15] assignments_to_experiments: keep the period that starts an experiment A period's time overrides the model's initial time, so an experiment starts at its first period after preequilibration. The converter dropped that period when nothing changed there, either as redundant, since the model already holds the value before the switch, or because its condition only set symbols that the conversion replaced. The experiment then started at its first switch, and measurements before the switch were simulated from the wrong start. The period that starts each experiment now stays, without conditions if none are left. A switch therefore never becomes an experiment's first period, and a start value that equals the model's no longer needs to be set when a switch refers to symbols outside the parameter table. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 39 ++++++----- tests/v2/test_conversion.py | 128 +++++++++++++++++++++++++++++------- 2 files changed, 127 insertions(+), 40 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 55bbdb92..38bacc2e 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -752,7 +752,10 @@ def assignments_to_experiments(problem: Problem) -> Problem: .. note:: Measurements without an experiment ID are assigned to a new experiment - :data:`DEFAULT_EXPERIMENT_ID`, if needed. + :data:`DEFAULT_EXPERIMENT_ID`, if needed, which starts at time 0. + Every experiment keeps the period that starts it, without conditions + if all of them were replaced, since a period's time overrides the + model's initial time. **Scope:** conversion is supported for SBML assignment rules that depend on time through the conditions of piecewise functions, and at most @@ -816,13 +819,22 @@ def assignments_to_experiments(problem: Problem) -> Problem: else: default_experiment = None # periods that are dropped if all their conditions were replaced - droppable = { - id(p) for e in experiments for p in e.periods if p.condition_ids - } + droppable = set() for experiment in experiments: if all(p.is_preequilibration for p in experiment.periods): - experiment += (period := ExperimentPeriod(time=0)) - droppable.add(id(period)) + experiment += ExperimentPeriod(time=0) + # The period that starts the experiment stays, even without + # conditions: its time overrides the model's initial time, so the + # experiment would start at its first switch instead. + start = min( + (p for p in experiment.periods if not p.is_preequilibration), + key=lambda p: p.time, + ) + droppable |= { + id(p) + for p in experiment.periods + if p.condition_ids and p is not start + } # (experiment, period start time, target ID, target value) additions = [] @@ -950,21 +962,16 @@ def assignments_to_experiments(problem: Problem) -> Problem: switch_values.append((t, value)) previous = value # The model already has the pre-switch value at the start of - # the first period, unless that value refers to parameter - # table symbols, for which the model only has a number. The - # switches can become the first periods if they only refer - # to parameter table symbols. A ramp's start and slope are - # always set by the tables. + # the first period, which stays even without conditions, + # unless that value refers to parameter table symbols, for + # which the model only has a number. A ramp's start and slope + # are always set by the tables. redundant = ( not ramps and i_period == 0 and (start_value - default).is_zero and not {str(s) for s in default.free_symbols} & parameter_ids - and all( - {str(s) for s in value.free_symbols} <= parameter_ids - for _, value in switch_values - ) ) after_preequilibration = ( i_period > 0 and periods[i_period - 1].is_preequilibration @@ -1049,7 +1056,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: p for p in experiment.periods if p.condition_ids or id(p) not in droppable - ] or experiment.periods[:1] + ] if default_experiment and any( p.condition_ids for p in default_experiment.periods diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index 306e63f9..f8b7f04a 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -14,6 +14,7 @@ from petab.v2.calculate import calculate_llh as calculate_llh_v2 from petab.v2.converters import ( DEFAULT_EXPERIMENT_ID, + ExperimentsToSbmlConverter, assignments_to_experiments, ) from petab.v2.petab1to2 import petab1to2, v1v2_observable_df @@ -685,14 +686,15 @@ def test_assignments_to_experiments(condition, before, after): converted = assignments_to_experiments(problem) # the value before the switch is kept in the model, so that only the - # switch itself needs a condition + # switch itself needs a condition; the experiment still starts at 0, + # where the first measurement is assert converted.model.sbml_model.getParameter("k1").getValue() == before experiments_expected = pd.DataFrame( data={ - v2.C.EXPERIMENT_ID: [DEFAULT_EXPERIMENT_ID], - v2.C.TIME: [10.0], - v2.C.CONDITION_ID: ["cond_0"], + v2.C.EXPERIMENT_ID: [DEFAULT_EXPERIMENT_ID] * 2, + v2.C.TIME: [0.0, 10.0], + v2.C.CONDITION_ID: ["", "cond_0"], } ) @@ -757,21 +759,22 @@ def test_assignments_to_experiments_nested(): assert converted.model.sbml_model.getParameter("k1").getValue() == 0.0 # `level` is not a parameter table entry, so it must not appear in the - # first period -- the pre-dose value is kept as a condition here + # first period -- the switch may refer to it, since the experiment + # starts with an empty period at 0 experiments_expected = pd.DataFrame( data={ v2.C.EXPERIMENT_ID: [DEFAULT_EXPERIMENT_ID] * 2, v2.C.TIME: [0.0, 10.0], - v2.C.CONDITION_ID: ["cond_0", "cond_1"], + v2.C.CONDITION_ID: ["", "cond_0"], } ) # the factor outside the piecewise expression is preserved conditions_expected = pd.DataFrame( data={ - v2.C.CONDITION_ID: ["cond_0", "cond_1"], - v2.C.TARGET_ID: ["k1", "k1"], - v2.C.TARGET_VALUE: [0.0, "level"], + v2.C.CONDITION_ID: ["cond_0"], + v2.C.TARGET_ID: ["k1"], + v2.C.TARGET_VALUE: ["level"], } ) @@ -865,14 +868,16 @@ def test_assignments_to_experiments_condition_specific(): converted = assignments_to_experiments(problem) # the pre-dose value is the same for both conditions and is kept in the - # model; each condition gets its own experiment and switching time + # model; each condition gets its own experiment and switching time, and + # each experiment still starts at 0, where its condition is now empty assert converted.model.sbml_model.getParameter("k1").getValue() == 0.0 experiments_expected = pd.DataFrame( data={ - v2.C.EXPERIMENT_ID: ["experiment____c1", "experiment____c2"], - v2.C.TIME: [10.0, 20.0], - v2.C.CONDITION_ID: ["cond_0", "cond_1"], + v2.C.EXPERIMENT_ID: ["experiment____c1"] * 2 + + ["experiment____c2"] * 2, + v2.C.TIME: [0.0, 10.0, 0.0, 20.0], + v2.C.CONDITION_ID: ["", "cond_0", "", "cond_1"], } ) @@ -893,6 +898,76 @@ def test_assignments_to_experiments_condition_specific(): ] +def _simulate(problem, experiment_id, timepoints): + """Simulate an experiment in BasiCO, after converting the experiments + to SBML events, and return the trajectory of `A`.""" + import basico + + converted = ExperimentsToSbmlConverter(problem).convert() + indicator = converted.model.sbml_model.getParameter( + ExperimentsToSbmlConverter.get_experiment_indicator(experiment_id) + ) + if indicator is not None: + indicator.setValue(1) + basico.load_model(converted.model.to_sbml_str()) + result = basico.run_time_course( + values=timepoints, r_tol=1e-10, a_tol=1e-12 + ) + return result["A"].to_numpy() + + +@pytest.mark.parametrize("by_condition", [False, True]) +def test_assignments_to_experiments_simulate_before_switch(by_condition): + """Test that a converted experiment starts at 0, not at its first switch, + if nothing changes at 0: without an experiment, or if the condition at + 0 only sets symbols of the converted assignment. Simulations before and + after the switch match the original problem's.""" + problem = Problem() + ant_model = f""" + model conversion + species A, B; + A = 10; + B = 0; + t_switch = 0; + k1 := piecewise(1, time < {"t_switch" if by_condition else 5}, 2); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + if by_condition: + experiment_id = "experiment____c1" + problem.add_condition("c1", t_switch=3.0) + problem.add_experiment(experiment_id, 0, "c1") + else: + experiment_id = DEFAULT_EXPERIMENT_ID + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + timepoints = [0.0, 1.0, 2.0, 4.0, 6.0, 10.0] + for t in timepoints: + problem.add_measurement( + "obs_a", + experiment_id=experiment_id if by_condition else None, + time=t, + measurement=5.0, + ) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + converted = assignments_to_experiments(problem) + + np.testing.assert_allclose( + _simulate(converted, experiment_id, timepoints), + _simulate(problem, experiment_id, timepoints), + rtol=1e-6, + ) + (experiment,) = converted.experiments + assert experiment.id == experiment_id + assert [(p.time, p.condition_ids) for p in experiment.periods] == [ + (0.0, []), + (3.0 if by_condition else 5.0, ["cond_0"]), + ] + + def test_assignments_to_experiments_preequilibration(): """Test that the converted assignments are applied during preequilibration, using the preequilibration condition's values.""" @@ -1040,13 +1115,16 @@ def test_assignments_to_experiments_fixed_switching_time( converted = assignments_to_experiments(problem) + periods = [ + (row[v2.C.TIME], row[v2.C.CONDITION_ID]) + for _, row in converted.experiment_df.iterrows() + ] + # the experiment starts at 0, where nothing changes + assert periods[0] == (0.0, "") conditions = converted.condition_df.set_index(v2.C.CONDITION_ID) assert [ - ( - row[v2.C.TIME], - float(conditions.loc[row[v2.C.CONDITION_ID], v2.C.TARGET_VALUE]), - ) - for _, row in converted.experiment_df.iterrows() + (t, float(conditions.loc[condition_id, v2.C.TARGET_VALUE])) + for t, condition_id in periods[1:] ] == [(t_switch, 1.0)] @@ -1136,14 +1214,16 @@ def test_assignments_to_experiments_several_switches( converted = assignments_to_experiments(problem) assert converted.model.sbml_model.getParameter("k1").getValue() == default - experiments = converted.experiment_df + periods = [ + (row[v2.C.TIME], row[v2.C.CONDITION_ID]) + for _, row in converted.experiment_df.iterrows() + ] + # the experiment starts at 0, where nothing changes + assert periods[0] == (0.0, "") conditions = converted.condition_df.set_index(v2.C.CONDITION_ID) assert [ - ( - row[v2.C.TIME], - float(conditions.loc[row[v2.C.CONDITION_ID], v2.C.TARGET_VALUE]), - ) - for _, row in experiments.iterrows() + (t, float(conditions.loc[condition_id, v2.C.TARGET_VALUE])) + for t, condition_id in periods[1:] ] == switches From 7bc856c15c0f16c14903402d0475152b2688f91e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 13:01:29 +0000 Subject: [PATCH 06/15] test_assignments_to_experiments_unused_condition_id: expect the start period Since the experiment keeps the period that starts it, the switch's condition is at the second period, after an empty one at t = 0. Co-Authored-By: Claude Opus 5.5 --- tests/v2/test_conversion.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index f8b7f04a..201a6740 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -729,7 +729,10 @@ def test_assignments_to_experiments_unused_condition_id(): converted = assignments_to_experiments(problem) assert [c.id for c in converted.conditions] == ["cond_0_1"] - assert converted.experiments[0].periods[0].condition_ids == ["cond_0_1"] + # the experiment starts at 0, where nothing changes + assert [ + (p.time, p.condition_ids) for p in converted.experiments[0].periods + ] == [(0.0, []), (10.0, ["cond_0_1"])] def test_assignments_to_experiments_nested(): From f8d2ac2a1051d0800bbc76b9c63a4583b2085b09 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Fri, 16 Oct 2026 09:10:16 +0100 Subject: [PATCH 07/15] assignments_to_experiments: convert values that ramp with time, via a rate rule A condition change is applied once at the start of its period, with `time` set to the period's start, so a value that changes with time within a period cannot be a condition change on a constant target, and such assignments were rejected. Fujita_SciSignal2010's EGF input is one: `piecewise(EGF_rate * time + EGF_0, time <= EGF_end, 0)`, steps in training and ramps in its held-out data. A value may now change linearly with time. The target then becomes a differential entity, with the rate rule `d target/dt = _slope` and a new constant parameter for the slope, and each period sets both the target's value at the period start and its slope. During preequilibration the value is held at its value at time 0 with slope zero, and the first period after preequilibration sets the slope again. An assignment that only switches between constant values converts exactly as before, as a constant parameter. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01MNJePorYs9zmPYMam6g4EP --- petab/v2/converters.py | 81 +++++++++++++++++--- tests/v2/test_conversion.py | 145 +++++++++++++++++++++++++++++++++++- 2 files changed, 213 insertions(+), 13 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 3c28a139..48ee2cdb 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -755,7 +755,8 @@ def assignments_to_experiments(problem: Problem) -> Problem: :data:`DEFAULT_EXPERIMENT_ID`, if needed. **Scope:** conversion is supported for SBML assignment rules that depend - on time only through the conditions of piecewise functions, where + on time through the conditions of piecewise functions, and at most + linearly through their values, where * every such condition is an inequality that can be re-arranged explicitly into `t />= ...`, @@ -770,7 +771,13 @@ def assignments_to_experiments(problem: Problem) -> Problem: given at different times (`d1 * piecewise(0, t < t1, 1) + d2 * piecewise(0, t < t2, 1)`), phases one after another as several pieces or nested piecewise functions (`piecewise(b0, t < t1, piecewise(b1, t < t2, - b2))`). Each switch that changes the value becomes a period. Anything else + b2))`). Each switch that changes the value becomes a period. A value may + also change linearly with time within a period (a ramp: `piecewise(r * t + + a, t <= t_end, 0)`). A condition change is applied once at the start of + its period, so the target of such an assignment becomes a differential + entity instead, with the rate rule `d target/dt = _slope`, and + each period sets both the target and its slope; an assignment that only + switches between constant values stays a constant parameter. Anything else raises a :class:`NotImplementedError`. .. warning:: @@ -855,7 +862,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: ) model_values = { s: sbml_model.getParameter(str(s)).getValue() - for s in default.free_symbols + for s in default.free_symbols - {time} } # unless a condition sets them, the switching times' symbols keep # their value from the parameter table or else from the model @@ -867,7 +874,24 @@ def assignments_to_experiments(problem: Problem) -> Problem: for s in t_switch.free_symbols } target = sbml_model.getParameter(target_id) - target.setValue(float(default.subs(model_values))) + # A value that changes with time within a period (a ramp) cannot be a + # condition change on a constant target, which is applied once at + # the start of its period. The target becomes a differential entity + # instead, `d target/dt = _slope`, and each period sets both. + ramps = expr.xreplace({c: sp.Dummy() for c in switches}).has(time) + slope_id = None + if ramps: + slope_id = f"{target_id}_slope" + while sbml_model.getElementBySId(slope_id) is not None: + slope_id += "_" + slope = sbml_model.createParameter() + slope.setId(slope_id) + slope.setConstant(False) + slope.setValue(float(sp.diff(default, time).subs(model_values))) + rate_rule = sbml_model.createRateRule() + rate_rule.setVariable(target_id) + rate_rule.setMath(libsbml.parseL3Formula(slope_id)) + target.setValue(float(default.subs(time, 0).subs(model_values))) target.setConstant(False) assignment.removeFromParentAndDelete() symbol_ids = {str(s) for s in expr.free_symbols} @@ -899,11 +923,16 @@ def assignments_to_experiments(problem: Problem) -> Problem: if period.is_preequilibration: # the value at time 0, which respects whether a - # condition holds for a switch at time 0 + # condition holds for a switch at time 0, held there: + # a ramp does not run during preequilibration value = expr.subs(subs).subs(time, 0) additions.append( (experiment, period.time, target_id, value) ) + if ramps: + additions.append( + (experiment, period.time, slope_id, sp.S.Zero) + ) continue t_end = periods[i + 1].time if i + 1 < len(periods) else sp.oo @@ -933,14 +962,27 @@ def assignments_to_experiments(problem: Problem) -> Problem: for _, value in switch_values ) ) + after_preequilibration = ( + i > 0 and periods[i - 1].is_preequilibration + ) if not redundant and ( - i == 0 or {c.target_id for c in changes} & symbol_ids + i == 0 + or {c.target_id for c in changes} & symbol_ids + # preequilibration held the slope at zero + or (ramps and after_preequilibration) ): - additions.append( - (experiment, period.time, target_id, start_value) + additions += _changes( + experiment, + period.time, + target_id, + start_value, + time, + slope_id, ) for t, value in switch_values: - additions.append((experiment, t, target_id, value)) + additions += _changes( + experiment, t, target_id, value, time, slope_id + ) if converted_ids and any(e.has_preequilibration for e in experiments): warnings.warn( @@ -1012,6 +1054,16 @@ def assignments_to_experiments(problem: Problem) -> Problem: return problem +def _changes(experiment, t, target_id, value, time, slope_id): + """The changes that start ``value`` at time ``t``: the target set to the + value at ``t``, and, for a target that ramps, its slope set to the + value's rate of change.""" + changes = [(experiment, t, target_id, value.subs(time, t))] + if slope_id is not None: + changes.append((experiment, t, slope_id, sp.diff(value, time))) + return changes + + def _value_from( expr: sp.Expr, switches: dict, @@ -1072,9 +1124,14 @@ def _supported_time_switches( raise unsupported t_switch = sp.solve(condition.gts - condition.lts, time)[0] switches[condition] = (t_switch, slope) - # time enters only through the conditions - if expr.xreplace({c: sp.Dummy() for c in switches}).has(time): - raise unsupported + # outside the conditions, time enters at most linearly: the value may + # ramp, at a rate that does not itself change with time + rest = expr.xreplace({c: sp.Dummy() for c in switches}) + if rest.has(time) and sp.diff(rest, time).has(time): + raise NotImplementedError( + "Only values that change at most linearly with time are " + f"supported, but the assignment to `{target_id}` is `{expr}`." + ) # no period can start at an estimated time if estimated := sorted( { diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index 7fd8fd96..ec9e7983 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -1,6 +1,7 @@ import logging from pathlib import Path +import libsbml import numpy as np import pandas as pd import pytest @@ -957,7 +958,8 @@ def test_assignments_to_experiments_preequilibration(): ("piecewise(0, time == 10, 1)", "single"), ("piecewise(0, time < 10 && time > 5, 1)", "single"), ("piecewise(0, time^2 < 10, 1)", "single"), - ("piecewise(time, time < 10, 1)", "single"), + ("piecewise(time^2, time < 10, 1)", "linearly"), + ("piecewise(exp(time), time < 10, 1)", "linearly"), ( "piecewise(0, time < t_dose, 1)", r"estimated parameters \['t_dose'\]", @@ -1196,6 +1198,147 @@ def test_assignments_to_experiments_condition_specific_doses(): } +def _rate_rule(problem, target_id): + rule = problem.model.sbml_model.getRateRuleByVariable(target_id) + return None if rule is None else libsbml.formulaToL3String(rule.getMath()) + + +def _changes_by_period(problem): + """{(experiment ID, time): {target ID: target value}}""" + conditions = problem.condition_df.set_index(v2.C.CONDITION_ID) + out = {} + for _, row in problem.experiment_df.iterrows(): + change = conditions.loc[[row[v2.C.CONDITION_ID]]] + out.setdefault((row[v2.C.EXPERIMENT_ID], row[v2.C.TIME]), {}).update( + { + r[v2.C.TARGET_ID]: float(r[v2.C.TARGET_VALUE]) + for _, r in change.iterrows() + } + ) + return out + + +def test_assignments_to_experiments_ramp(): + """Test that a value changing linearly with time becomes a rate rule + whose slope each period sets.""" + problem = Problem() + ant_model = """ + model conversion + species A, B; + A = 10; + B = 0; + r = 0.5; a = 1; + k1 := piecewise(r * time + a, time <= 10, 0); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement("obs_a", time=20.0, measurement=5.0) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + converted = assignments_to_experiments(problem) + + sbml = converted.model.sbml_model + # the ramp's start and slope are the model's own initial values + assert _rate_rule(converted, "k1") == "k1_slope" + assert sbml.getParameter("k1").getValue() == 1.0 + assert sbml.getParameter("k1_slope").getValue() == 0.5 + # at the switch, the target and its slope drop to zero + assert _changes_by_period(converted) == { + (DEFAULT_EXPERIMENT_ID, 10.0): {"k1": 0.0, "k1_slope": 0.0} + } + + +def test_assignments_to_experiments_condition_specific_ramps(): + """Test pulses and ramps whose start, rate and end are set by the + condition table, as in Fujita_SciSignal2010.""" + problem = Problem() + ant_model = """ + model conversion + species A, B; + A = 10; + B = 0; + EGF_0 = 0; EGF_rate = 0; EGF_end = 0; + k1 := piecewise(EGF_rate * time + EGF_0, time <= EGF_end, 0); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + # a pulse: a constant dose withdrawn at t = 60 + problem.add_condition("pulse", EGF_0=1.0, EGF_rate=0.0, EGF_end=60.0) + # a ramp: growing from zero at 0.01 per unit time until t = 3600 + problem.add_condition("ramp", EGF_0=0.0, EGF_rate=0.01, EGF_end=3600.0) + for c in ("pulse", "ramp"): + problem.add_experiment(f"experiment____{c}", 0, c) + problem.add_measurement( + "obs_a", + experiment_id=f"experiment____{c}", + time=30.0, + measurement=5.0, + ) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + converted = assignments_to_experiments(problem) + + assert _rate_rule(converted, "k1") == "k1_slope" + assert _changes_by_period(converted) == { + ("experiment____pulse", 0.0): {"k1": 1.0, "k1_slope": 0.0}, + ("experiment____pulse", 60.0): {"k1": 0.0, "k1_slope": 0.0}, + ("experiment____ramp", 0.0): {"k1": 0.0, "k1_slope": 0.01}, + ("experiment____ramp", 3600.0): {"k1": 0.0, "k1_slope": 0.0}, + } + + +def test_assignments_to_experiments_ramp_preequilibration(): + """Test that a ramp is held at its value at time 0 during + preequilibration and starts at time 0.""" + problem = Problem() + ant_model = """ + model conversion + species A, B; + A = 10; + B = 0; + r = 0; + k1 := piecewise(r * time + 1, time <= 10, 0); + k2 = 0.5; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; + end + """ + problem.model = SbmlModel.from_antimony(ant_model) + problem.add_condition("pre", r=0.0) + problem.add_condition("main", r=2.0) + problem.add_experiment( + "experiment__pre___main", v2.C.TIME_PREEQUILIBRATION, "pre", 0, "main" + ) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", + experiment_id="experiment__pre___main", + time=20.0, + measurement=5.0, + ) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + + with pytest.warns(UserWarning, match="during preequilibration"): + converted = assignments_to_experiments(problem) + + assert _changes_by_period(converted) == { + ("experiment__pre___main", v2.C.TIME_PREEQUILIBRATION): { + "k1": 1.0, + "k1_slope": 0.0, + }, + ("experiment__pre___main", 0.0): {"k1": 1.0, "k1_slope": 2.0}, + ("experiment__pre___main", 10.0): {"k1": 0.0, "k1_slope": 0.0}, + } + + @pytest.mark.parametrize("t_switch", [0, 10]) @pytest.mark.parametrize( "condition, holds_at_switch, holds_before", From 86d35a92446344a067f947158ae326dc8bb72f28 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Fri, 16 Oct 2026 09:49:12 +0100 Subject: [PATCH 08/15] assignments_to_experiments: model change only makes a ramp's target differential The slope parameter is a placeholder at 0; every ramp's start value and slope, including at time 0, are set through the condition table. Constant model parameters in changes at an experiment's start take their model values, since those changes may only refer to parameter table symbols. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_01MNJePorYs9zmPYMam6g4EP --- petab/v2/converters.py | 24 +++++++++++++++++++----- tests/v2/test_conversion.py | 12 +++++++----- 2 files changed, 26 insertions(+), 10 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 48ee2cdb..0ef9a14f 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -876,8 +876,9 @@ def assignments_to_experiments(problem: Problem) -> Problem: target = sbml_model.getParameter(target_id) # A value that changes with time within a period (a ramp) cannot be a # condition change on a constant target, which is applied once at - # the start of its period. The target becomes a differential entity - # instead, `d target/dt = _slope`, and each period sets both. + # the start of its period. The model only makes the target a + # differential entity, `d target/dt = _slope`; the slope's + # value, like the target's, is set by every period that starts one. ramps = expr.xreplace({c: sp.Dummy() for c in switches}).has(time) slope_id = None if ramps: @@ -887,7 +888,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: slope = sbml_model.createParameter() slope.setId(slope_id) slope.setConstant(False) - slope.setValue(float(sp.diff(default, time).subs(model_values))) + slope.setValue(0.0) rate_rule = sbml_model.createRateRule() rate_rule.setVariable(target_id) rate_rule.setMath(libsbml.parseL3Formula(slope_id)) @@ -951,9 +952,11 @@ def assignments_to_experiments(problem: Problem) -> Problem: # the first period, unless that value refers to parameter # table symbols, for which the model only has a number. The # switches can become the first periods if they only refer - # to parameter table symbols. + # to parameter table symbols. A ramp's start and slope are + # always set by the tables. redundant = ( - i == 0 + not ramps + and i == 0 and (start_value - default).is_zero and not {str(s) for s in default.free_symbols} & parameter_ids @@ -971,6 +974,17 @@ def assignments_to_experiments(problem: Problem) -> Problem: # preequilibration held the slope at zero or (ramps and after_preequilibration) ): + if i == 0: + # at the start of an experiment, a change may only + # refer to parameter table symbols; other symbols + # are constant model parameters, at their values + start_value = start_value.subs( + { + s: sbml_model.getParameter(str(s)).getValue() + for s in start_value.free_symbols - {time} + if str(s) not in parameter_ids + } + ) additions += _changes( experiment, period.time, diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index ec9e7983..306e63f9 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -1242,13 +1242,15 @@ def test_assignments_to_experiments_ramp(): converted = assignments_to_experiments(problem) sbml = converted.model.sbml_model - # the ramp's start and slope are the model's own initial values + # the model only makes the target differential; the slope is a + # placeholder until the tables set it assert _rate_rule(converted, "k1") == "k1_slope" - assert sbml.getParameter("k1").getValue() == 1.0 - assert sbml.getParameter("k1_slope").getValue() == 0.5 - # at the switch, the target and its slope drop to zero + assert sbml.getParameter("k1_slope").getValue() == 0.0 + # the ramp starts at time 0, and at the switch, the target and its slope + # drop to zero assert _changes_by_period(converted) == { - (DEFAULT_EXPERIMENT_ID, 10.0): {"k1": 0.0, "k1_slope": 0.0} + (DEFAULT_EXPERIMENT_ID, 0.0): {"k1": 1.0, "k1_slope": 0.5}, + (DEFAULT_EXPERIMENT_ID, 10.0): {"k1": 0.0, "k1_slope": 0.0}, } From 1c90569485af743586e86d6762b62de15a3b3ccf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Mon, 19 Oct 2026 19:32:12 +0100 Subject: [PATCH 09/15] assignments_to_experiments: rename the period index to `i_period` Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 0ef9a14f..7dcd5a66 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -902,7 +902,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: for experiment in experiments: periods = experiment.sorted_periods overrides = {} - for i, period in enumerate(periods): + for i_period, period in enumerate(periods): changes = problem.get_changes_for_period(period) overrides |= {c.target_id: c.target_value for c in changes} subs = fixed_values | { @@ -936,7 +936,11 @@ def assignments_to_experiments(problem: Problem) -> Problem: ) continue - t_end = periods[i + 1].time if i + 1 < len(periods) else sp.oo + t_end = ( + periods[i_period + 1].time + if i_period + 1 < len(periods) + else sp.oo + ) start_value = _value_from( expr, switches, time, subs, times, period.time ) @@ -956,7 +960,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: # always set by the tables. redundant = ( not ramps - and i == 0 + and i_period == 0 and (start_value - default).is_zero and not {str(s) for s in default.free_symbols} & parameter_ids @@ -966,15 +970,15 @@ def assignments_to_experiments(problem: Problem) -> Problem: ) ) after_preequilibration = ( - i > 0 and periods[i - 1].is_preequilibration + i_period > 0 and periods[i_period - 1].is_preequilibration ) if not redundant and ( - i == 0 + i_period == 0 or {c.target_id for c in changes} & symbol_ids # preequilibration held the slope at zero or (ramps and after_preequilibration) ): - if i == 0: + if i_period == 0: # at the start of an experiment, a change may only # refer to parameter table symbols; other symbols # are constant model parameters, at their values From c52a56fae5935d005d60e1b73222ebe541018824 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 12:53:12 +0000 Subject: [PATCH 10/15] assignments_to_experiments: add a ramp's slope with `add_sbml_parameter` `add_sbml_parameter` checks libsbml's return codes, which the hand-written parameter creation did not. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 7dcd5a66..55bbdb92 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -885,10 +885,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: slope_id = f"{target_id}_slope" while sbml_model.getElementBySId(slope_id) is not None: slope_id += "_" - slope = sbml_model.createParameter() - slope.setId(slope_id) - slope.setConstant(False) - slope.setValue(0.0) + add_sbml_parameter(sbml_model, slope_id, value=0.0, constant=False) rate_rule = sbml_model.createRateRule() rate_rule.setVariable(target_id) rate_rule.setMath(libsbml.parseL3Formula(slope_id)) From 50cac56e845cfc151b8aaa8cf083f4f4946609fb Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Mon, 19 Oct 2026 19:38:15 +0100 Subject: [PATCH 11/15] assignments_to_experiments: split into smaller pieces The function had grown to 340 lines, most of them one loop body passing around the parsed assignment, its switches, its value before every switch, the switching times' fixed values and the slope parameter of a ramp. These now live on `_PiecewiseAssignment`: its constructor checks that the assignment is supported, `replace_in_model` replaces it in the model, and `changes` returns the condition changes that replace it in an experiment. `_changes` and `_value_from` become its methods. The "redundant" check becomes a decision per period, one branch for the first period and one for the later ones. Setting up the start periods, dropping the replaced changes and adding the new conditions are separate functions. No change in behaviour: the old and new implementations give the same tables, SBML, errors and warnings on all conversions in the tests, and on the benchmark problems Brannmark_JBC2010, Fujita_SciSignal2010, Isensee_JCB2018, Oliveira_NatCommun2021 and Raimundez_PCB2020. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 502 +++++++++++++++++++++++------------------ 1 file changed, 277 insertions(+), 225 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 38bacc2e..9eef1061 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -806,11 +806,6 @@ def assignments_to_experiments(problem: Problem) -> Problem: problem = deepcopy(problem) sbml_model = problem.model.sbml_model - parameter_ids = {p.id for p in problem.parameters} - estimated_ids = {p.id for p in problem.parameters if p.estimate} - nominal_values = { - p.id: p.nominal_value for p in problem.parameters if not p.estimate - } experiments = list(problem.experiments) if any(m.experiment_id is None for m in problem.measurements): @@ -818,224 +813,332 @@ def assignments_to_experiments(problem: Problem) -> Problem: experiments.append(default_experiment) else: default_experiment = None - # periods that are dropped if all their conditions were replaced - droppable = set() - for experiment in experiments: - if all(p.is_preequilibration for p in experiment.periods): - experiment += ExperimentPeriod(time=0) - # The period that starts the experiment stays, even without - # conditions: its time overrides the model's initial time, so the - # experiment would start at its first switch instead. - start = min( - (p for p in experiment.periods if not p.is_preequilibration), - key=lambda p: p.time, - ) - droppable |= { - id(p) - for p in experiment.periods - if p.condition_ids and p is not start - } + droppable = _add_start_periods(experiments) # (experiment, period start time, target ID, target value) additions = [] replaced_symbols = set() converted_ids = [] - - for assignment in [ + for rule in [ *(r for r in sbml_model.getListOfRules() if r.isAssignment()), *sbml_model.getListOfInitialAssignments(), ]: - expr = sbml_math_to_sympy(assignment) + expr = sbml_math_to_sympy(rule) if not any(pw.has(TimeSymbol) for pw in expr.atoms(sp.Piecewise)): continue - target_id = assignment.getId() - (time,) = expr.atoms(TimeSymbol) - switches = _supported_time_switches( - expr, time, target_id, estimated_ids + assignment = _PiecewiseAssignment(rule, expr, problem) + assignment.replace_in_model() + for experiment in experiments: + additions += [ + (experiment, *change) + for change in assignment.changes(experiment) + ] + replaced_symbols |= assignment.symbol_ids + converted_ids.append(assignment.target_id) + + if converted_ids and any(e.has_preequilibration for e in experiments): + warnings.warn( + "Switches of time-dependent piecewise assignments during " + "preequilibration are not supported. The converted " + f"assignments to {converted_ids} keep their value at time 0 " + "throughout " + "preequilibration.", + stacklevel=2, + ) + + # replace the changes that are covered by the new conditions + _drop_changes(problem, experiments, replaced_symbols) + _add_conditions(problem, additions) + + for experiment in experiments: + experiment.sort_periods() + experiment.periods = [ + p + for p in experiment.periods + if p.condition_ids or id(p) not in droppable + ] + + if default_experiment and any( + p.condition_ids for p in default_experiment.periods + ): + problem += default_experiment + for measurement in problem.measurements: + if measurement.experiment_id is None: + measurement.experiment_id = DEFAULT_EXPERIMENT_ID + + validation_issues = problem.validate() + if validation_issues.has_errors(): + messages = "\n".join(f"{issue.message}" for issue in validation_issues) + raise ValueError( + f"The converted PEtab problem did not pass linting:\n{messages}" + ) + + return problem + + +class _PiecewiseAssignment: + """A time-dependent piecewise assignment to a model parameter, and the + condition changes that replace it, see :func:`assignments_to_experiments`. + """ + + def __init__( + self, + rule: libsbml.Rule | libsbml.InitialAssignment, + expr: sp.Expr, + problem: Problem, + ): + """Check that the assignment ``rule``, with the value ``expr``, is + supported.""" + self.rule = rule + self.problem = problem + self.sbml_model = problem.model.sbml_model + self.target_id = rule.getId() + self.expr = expr + (self.time,) = expr.atoms(TimeSymbol) + self.symbol_ids = {str(s) for s in expr.free_symbols} + self.parameter_ids = {p.id for p in problem.parameters} + self.switches = _supported_time_switches( + expr, + self.time, + self.target_id, + {p.id for p in problem.parameters if p.estimate}, ) if non_constant := sorted( str(s) - for s in expr.free_symbols - {time} - if not (p := sbml_model.getParameter(str(s))) + for s in expr.free_symbols - {self.time} + if not (p := self.sbml_model.getParameter(str(s))) or not p.getConstant() ): raise NotImplementedError( "Only constant model parameters are supported as symbols in " - f"the assignment to `{target_id}`, but it refers to " + f"the assignment to `{self.target_id}`, but it refers to " f"{non_constant}." ) # the value before every switch: a condition whose two sides grow # apart with time holds only after its switch - default = expr.xreplace( + self.default = expr.xreplace( { c: sp.false if slope > 0 else sp.true - for c, (_, slope) in switches.items() + for c, (_, slope) in self.switches.items() } ) - model_values = { - s: sbml_model.getParameter(str(s)).getValue() - for s in default.free_symbols - {time} - } # unless a condition sets them, the switching times' symbols keep # their value from the parameter table or else from the model - fixed_values = { - s: nominal_values.get( - str(s), sbml_model.getParameter(str(s)).getValue() - ) - for t_switch, _ in switches.values() + nominal_values = { + p.id: p.nominal_value for p in problem.parameters if not p.estimate + } + self.fixed_values = { + s: nominal_values.get(str(s), self._model_value(s)) + for t_switch, _ in self.switches.values() for s in t_switch.free_symbols } - target = sbml_model.getParameter(target_id) - # A value that changes with time within a period (a ramp) cannot be a - # condition change on a constant target, which is applied once at - # the start of its period. The model only makes the target a - # differential entity, `d target/dt = _slope`; the slope's - # value, like the target's, is set by every period that starts one. - ramps = expr.xreplace({c: sp.Dummy() for c in switches}).has(time) - slope_id = None - if ramps: - slope_id = f"{target_id}_slope" - while sbml_model.getElementBySId(slope_id) is not None: - slope_id += "_" - add_sbml_parameter(sbml_model, slope_id, value=0.0, constant=False) - rate_rule = sbml_model.createRateRule() - rate_rule.setVariable(target_id) - rate_rule.setMath(libsbml.parseL3Formula(slope_id)) - target.setValue(float(default.subs(time, 0).subs(model_values))) + # whether the value changes with time within a period + self.ramps = expr.xreplace({c: sp.Dummy() for c in self.switches}).has( + self.time + ) + self.slope_id = None + + def replace_in_model(self) -> None: + """Replace the assignment in the model by a target that the + condition table sets. + + A value that changes with time within a period (a ramp) cannot be a + condition change on a constant target, which is applied once at the + start of its period. The model only makes the target a differential + entity, `d target/dt = _slope`; the slope's value, like the + target's, is set by every period that starts one. + """ + if self.ramps: + self.slope_id = f"{self.target_id}_slope" + while self.sbml_model.getElementBySId(self.slope_id) is not None: + self.slope_id += "_" + add_sbml_parameter( + self.sbml_model, self.slope_id, value=0.0, constant=False + ) + rate_rule = self.sbml_model.createRateRule() + rate_rule.setVariable(self.target_id) + rate_rule.setMath(libsbml.parseL3Formula(self.slope_id)) + model_values = { + s: self._model_value(s) + for s in self.default.free_symbols - {self.time} + } + target = self.sbml_model.getParameter(self.target_id) + target.setValue( + float(self.default.subs(self.time, 0).subs(model_values)) + ) target.setConstant(False) - assignment.removeFromParentAndDelete() - symbol_ids = {str(s) for s in expr.free_symbols} - replaced_symbols |= symbol_ids - converted_ids.append(target_id) + self.rule.removeFromParentAndDelete() + + def changes( + self, experiment: Experiment + ) -> list[tuple[float, str, sp.Expr]]: + """The changes that replace the assignment in ``experiment``, each as + ``(period start time, target ID, target value)``.""" + changes = [] + periods = experiment.sorted_periods + overrides = {} + for i_period, period in enumerate(periods): + period_changes = self.problem.get_changes_for_period(period) + overrides |= {c.target_id: c.target_value for c in period_changes} + subs = self.fixed_values | { + s: overrides[str(s)] + for s in self.expr.free_symbols + if str(s) in overrides + } + times = self._switching_times(subs, experiment.id) - for experiment in experiments: - periods = experiment.sorted_periods - overrides = {} - for i_period, period in enumerate(periods): - changes = problem.get_changes_for_period(period) - overrides |= {c.target_id: c.target_value for c in changes} - subs = fixed_values | { - s: overrides[str(s)] - for s in expr.free_symbols - if str(s) in overrides - } - times = [] - for t_switch, _ in switches.values(): - t = t_switch.subs(subs) - if not (t.is_number and t.is_real and t.is_finite): - raise NotImplementedError( - "Only switching times that evaluate to a real " - "number are supported, but a switching time of " - f"the assignment to `{target_id}` in experiment " - f"`{experiment.id}` is `{t}`." - ) - times.append(float(t)) - - if period.is_preequilibration: - # the value at time 0, which respects whether a - # condition holds for a switch at time 0, held there: - # a ramp does not run during preequilibration - value = expr.subs(subs).subs(time, 0) - additions.append( - (experiment, period.time, target_id, value) - ) - if ramps: - additions.append( - (experiment, period.time, slope_id, sp.S.Zero) - ) - continue - - t_end = ( - periods[i_period + 1].time - if i_period + 1 < len(periods) - else sp.oo - ) - start_value = _value_from( - expr, switches, time, subs, times, period.time - ) - # each switch within the period that changes the value - switch_values = [] - previous = start_value - for t in sorted({t for t in times if period.time < t < t_end}): - value = _value_from(expr, switches, time, subs, times, t) - if not (value - previous).is_zero: - switch_values.append((t, value)) - previous = value + if period.is_preequilibration: + # the value at time 0, which respects whether a condition + # holds for a switch at time 0, held there: a ramp does not + # run during preequilibration + value = self.expr.subs(subs).subs(self.time, 0) + changes.append((period.time, self.target_id, value)) + if self.ramps: + changes.append((period.time, self.slope_id, sp.S.Zero)) + continue + + start_value = self._value_from(subs, times, period.time) + # each switch within the period that changes the value + t_end = ( + periods[i_period + 1].time + if i_period + 1 < len(periods) + else sp.oo + ) + switch_values = [] + previous = start_value + for t in sorted({t for t in times if period.time < t < t_end}): + value = self._value_from(subs, times, t) + if not (value - previous).is_zero: + switch_values.append((t, value)) + previous = value + + if i_period == 0: # The model already has the pre-switch value at the start of # the first period, which stays even without conditions, # unless that value refers to parameter table symbols, for - # which the model only has a number. A ramp's start and slope - # are always set by the tables. - redundant = ( - not ramps - and i_period == 0 - and (start_value - default).is_zero - and not {str(s) for s in default.free_symbols} - & parameter_ids - ) - after_preequilibration = ( - i_period > 0 and periods[i_period - 1].is_preequilibration - ) - if not redundant and ( - i_period == 0 - or {c.target_id for c in changes} & symbol_ids - # preequilibration held the slope at zero - or (ramps and after_preequilibration) + # which the model only has a number. A ramp's start and + # slope are always set by the tables. + if ( + self.ramps + or not (start_value - self.default).is_zero + or {str(s) for s in self.default.free_symbols} + & self.parameter_ids ): - if i_period == 0: - # at the start of an experiment, a change may only - # refer to parameter table symbols; other symbols - # are constant model parameters, at their values - start_value = start_value.subs( - { - s: sbml_model.getParameter(str(s)).getValue() - for s in start_value.free_symbols - {time} - if str(s) not in parameter_ids - } - ) - additions += _changes( - experiment, - period.time, - target_id, - start_value, - time, - slope_id, + # at the start of an experiment, a change may only refer + # to parameter table symbols; other symbols are + # constant model parameters, at their values + start_value = start_value.subs( + { + s: self._model_value(s) + for s in start_value.free_symbols - {self.time} + if str(s) not in self.parameter_ids + } ) - for t, value in switch_values: - additions += _changes( - experiment, t, target_id, value, time, slope_id - ) - - if converted_ids and any(e.has_preequilibration for e in experiments): - warnings.warn( - "Switches of time-dependent piecewise assignments during " - "preequilibration are not supported. The converted " - f"assignments to {converted_ids} keep their value at time 0 " - "throughout " - "preequilibration.", - stacklevel=2, + changes += self._start(period.time, start_value) + elif ( + {c.target_id for c in period_changes} & self.symbol_ids + # preequilibration held the slope at zero + or (self.ramps and periods[i_period - 1].is_preequilibration) + ): + changes += self._start(period.time, start_value) + for t, value in switch_values: + changes += self._start(t, value) + return changes + + def _model_value(self, symbol: sp.Symbol) -> float: + """The value of the model parameter ``symbol``.""" + return self.sbml_model.getParameter(str(symbol)).getValue() + + def _switching_times(self, subs: dict, experiment_id: str) -> list[float]: + """The switching times, for the symbol values ``subs``.""" + times = [] + for t_switch, _ in self.switches.values(): + t = t_switch.subs(subs) + if not (t.is_number and t.is_real and t.is_finite): + raise NotImplementedError( + "Only switching times that evaluate to a real number are " + "supported, but a switching time of the assignment to " + f"`{self.target_id}` in experiment `{experiment_id}` is " + f"`{t}`." + ) + times.append(float(t)) + return times + + def _value_from(self, subs: dict, times: list[float], t0: float): + """The value from ``t0`` up to its next switch, for the symbol values + ``subs``, every condition replaced by whether it holds in between.""" + later = [t for t in times if t > t0] + probe = (t0 + min(later)) / 2 if later else t0 + 1 + truth = {c: c.subs(subs).subs(self.time, probe) for c in self.switches} + return self.expr.xreplace(truth).subs(subs) + + def _start( + self, t: float, value: sp.Expr + ) -> list[tuple[float, str, sp.Expr]]: + """The changes that start ``value`` at time ``t``: the target set to + the value at ``t``, and, for a target that ramps, its slope set to + the value's rate of change.""" + changes = [(t, self.target_id, value.subs(self.time, t))] + if self.ramps: + changes.append((t, self.slope_id, sp.diff(value, self.time))) + return changes + + +def _add_start_periods(experiments: list[Experiment]) -> set[int]: + """Give every experiment a period that starts it, and return the + :func:`id` of the periods that are dropped if all their conditions are + replaced. + + The period that starts an experiment stays, even without conditions: its + time overrides the model's initial time, so the experiment would start at + its first switch instead. + """ + droppable = set() + for experiment in experiments: + if all(p.is_preequilibration for p in experiment.periods): + experiment += ExperimentPeriod(time=0) + start = min( + (p for p in experiment.periods if not p.is_preequilibration), + key=lambda p: p.time, ) + droppable |= { + id(p) + for p in experiment.periods + if p.condition_ids and p is not start + } + return droppable - # drop changes that are covered by the new conditions + +def _drop_changes( + problem: Problem, experiments: list[Experiment], target_ids: set[str] +) -> None: + """Drop the changes of ``target_ids`` from the condition tables, and the + conditions left without changes from the periods of ``experiments``.""" condition_ids = {c.id for c in problem.conditions} for table in problem.condition_tables: for condition in table.conditions: condition.changes = [ - c - for c in condition.changes - if c.target_id not in replaced_symbols + c for c in condition.changes if c.target_id not in target_ids ] table.elements = [c for c in table.conditions if c.changes] - removed_condition_ids = condition_ids - {c.id for c in problem.conditions} + removed_ids = condition_ids - {c.id for c in problem.conditions} for experiment in experiments: for period in experiment.periods: period.condition_ids = [ - cid - for cid in period.condition_ids - if cid not in removed_condition_ids + cid for cid in period.condition_ids if cid not in removed_ids ] - # add the new conditions to new or existing periods + +def _add_conditions( + problem: Problem, + additions: list[tuple[Experiment, float, str, sp.Expr]], +) -> None: + """For each ``(experiment, time, target ID, value)`` in ``additions``, + add a condition with that change to the period of ``experiment`` that + starts at ``time``, which is added if needed. + + Equal changes share a condition. + """ new_conditions = {} for experiment, t, target_id, value in additions: if (cond_id := new_conditions.get((target_id, value))) is None: @@ -1050,57 +1153,6 @@ def assignments_to_experiments(problem: Problem) -> Problem: experiment += (period := ExperimentPeriod(time=t)) period.condition_ids.append(cond_id) - for experiment in experiments: - experiment.sort_periods() - experiment.periods = [ - p - for p in experiment.periods - if p.condition_ids or id(p) not in droppable - ] - - if default_experiment and any( - p.condition_ids for p in default_experiment.periods - ): - problem += default_experiment - for measurement in problem.measurements: - if measurement.experiment_id is None: - measurement.experiment_id = DEFAULT_EXPERIMENT_ID - - validation_issues = problem.validate() - if validation_issues.has_errors(): - messages = "\n".join(f"{issue.message}" for issue in validation_issues) - raise ValueError( - f"The converted PEtab problem did not pass linting:\n{messages}" - ) - - return problem - - -def _changes(experiment, t, target_id, value, time, slope_id): - """The changes that start ``value`` at time ``t``: the target set to the - value at ``t``, and, for a target that ramps, its slope set to the - value's rate of change.""" - changes = [(experiment, t, target_id, value.subs(time, t))] - if slope_id is not None: - changes.append((experiment, t, slope_id, sp.diff(value, time))) - return changes - - -def _value_from( - expr: sp.Expr, - switches: dict, - time: sp.Symbol, - subs: dict, - times: list[float], - t0: float, -) -> sp.Expr: - """The value of ``expr`` from ``t0`` up to its next switch, every - condition replaced by whether it holds in between.""" - later = [t for t in times if t > t0] - probe = (t0 + min(later)) / 2 if later else t0 + 1 - truth = {c: c.subs(subs).subs(time, probe) for c in switches} - return expr.xreplace(truth).subs(subs) - def _supported_time_switches( expr: sp.Expr, time: sp.Symbol, target_id: str, estimated_ids: set[str] From 81d8742bcdf29ea2eb63c12c59758208c2efed5c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 12:26:49 +0000 Subject: [PATCH 12/15] assignments_to_experiments: fix lost switches, dropped changes and scope gaps - A switch exactly at the start of a later period was lost when that period did not change any symbol of the assignment, e.g. `piecewise(0, time < 10, 1)` with a period at t = 10, or `time <= 0` after preequilibration. A later period now also sets the value when it differs from the value at the end of the previous period. - The condition changes of every symbol of a converted assignment were dropped, even when the model still used the symbol elsewhere, e.g. a dose that also enters a reaction. They are now only dropped if nothing else (the model, an observable, a condition change) refers to the symbol. - Time-dependent initial assignments were converted into switches, but SBML evaluates them once, at the start of the simulation. They now stay in the model. - At preequilibration, which starts the experiment, a value that referred to model parameters outside the parameter table failed the lint. Those parameters now take their model values, as at the first period. - An assignment to a species or compartment crashed with an AttributeError; it now raises NotImplementedError. - A problem without time-dependent assignments was still changed: an experiment that only preequilibrates got an empty period at t = 0. `test_assignments_to_experiments` and `_unused_condition_id` used an initial assignment to mean a switching rule and now use an assignment rule. `_validate` relied on an estimated parameter with a time-dependent initial assignment; it now uses a ramp whose generated slope parameter collides with an observable parameter. The converted benchmark problems (Brannmark_JBC2010, Fujita_SciSignal2010, Isensee_JCB2018, Raimundez_PCB2020) are unchanged. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 115 ++++++++++++++++----- tests/v2/test_conversion.py | 195 ++++++++++++++++++++++++++++++++++-- 2 files changed, 278 insertions(+), 32 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index 9eef1061..b0f5b0b5 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -757,9 +757,9 @@ def assignments_to_experiments(problem: Problem) -> Problem: if all of them were replaced, since a period's time overrides the model's initial time. - **Scope:** conversion is supported for SBML assignment rules that depend - on time through the conditions of piecewise functions, and at most - linearly through their values, where + **Scope:** conversion is supported for SBML assignment rules to + parameters that depend on time through the conditions of piecewise + functions, and at most linearly through their values, where * every such condition is an inequality that can be re-arranged explicitly into `t />= ...`, @@ -781,7 +781,11 @@ def assignments_to_experiments(problem: Problem) -> Problem: entity instead, with the rate rule `d target/dt = _slope`, and each period sets both the target and its slope; an assignment that only switches between constant values stays a constant parameter. Anything else - raises a :class:`NotImplementedError`. + raises a :class:`NotImplementedError`. Initial assignments stay in the + model, since they are evaluated once, at the start of the simulation. + + Condition changes of the replaced assignments' symbols are dropped, unless + the model, an observable or another change still refers to them. .. warning:: @@ -807,25 +811,29 @@ def assignments_to_experiments(problem: Problem) -> Problem: problem = deepcopy(problem) sbml_model = problem.model.sbml_model + # Initial assignments are evaluated once, at the start of the + # simulation, so they never switch and stay in the model. + time_dependent = [] + for rule in sbml_model.getListOfRules(): + if not rule.isAssignment(): + continue + expr = sbml_math_to_sympy(rule) + if any(pw.has(TimeSymbol) for pw in expr.atoms(sp.Piecewise)): + time_dependent.append((rule, expr)) + experiments = list(problem.experiments) if any(m.experiment_id is None for m in problem.measurements): default_experiment = Experiment(id=DEFAULT_EXPERIMENT_ID) experiments.append(default_experiment) else: default_experiment = None - droppable = _add_start_periods(experiments) + droppable = _add_start_periods(experiments) if time_dependent else set() # (experiment, period start time, target ID, target value) additions = [] replaced_symbols = set() converted_ids = [] - for rule in [ - *(r for r in sbml_model.getListOfRules() if r.isAssignment()), - *sbml_model.getListOfInitialAssignments(), - ]: - expr = sbml_math_to_sympy(rule) - if not any(pw.has(TimeSymbol) for pw in expr.atoms(sp.Piecewise)): - continue + for rule, expr in time_dependent: assignment = _PiecewiseAssignment(rule, expr, problem) assignment.replace_in_model() for experiment in experiments: @@ -846,7 +854,9 @@ def assignments_to_experiments(problem: Problem) -> Problem: stacklevel=2, ) - # replace the changes that are covered by the new conditions + # replace the changes that are covered by the new conditions, unless + # something else still refers to their targets + replaced_symbols -= _referenced_ids(problem, additions) _drop_changes(problem, experiments, replaced_symbols) _add_conditions(problem, additions) @@ -893,6 +903,11 @@ def __init__( self.problem = problem self.sbml_model = problem.model.sbml_model self.target_id = rule.getId() + if self.sbml_model.getParameter(self.target_id) is None: + raise NotImplementedError( + "Only assignments to parameters are supported, but " + f"`{self.target_id}` is not a parameter." + ) self.expr = expr (self.time,) = expr.atoms(TimeSymbol) self.symbol_ids = {str(s) for s in expr.free_symbols} @@ -977,6 +992,8 @@ def changes( changes = [] periods = experiment.sorted_periods overrides = {} + # the value at the end of the previous period + held = None for i_period, period in enumerate(periods): period_changes = self.problem.get_changes_for_period(period) overrides |= {c.target_id: c.target_value for c in period_changes} @@ -991,8 +1008,14 @@ def changes( # the value at time 0, which respects whether a condition # holds for a switch at time 0, held there: a ramp does not # run during preequilibration - value = self.expr.subs(subs).subs(self.time, 0) - changes.append((period.time, self.target_id, value)) + held = self.expr.subs(subs).subs(self.time, 0) + changes.append( + ( + period.time, + self.target_id, + self._with_model_values(held), + ) + ) if self.ramps: changes.append((period.time, self.slope_id, sp.S.Zero)) continue @@ -1024,31 +1047,38 @@ def changes( or {str(s) for s in self.default.free_symbols} & self.parameter_ids ): - # at the start of an experiment, a change may only refer - # to parameter table symbols; other symbols are - # constant model parameters, at their values - start_value = start_value.subs( - { - s: self._model_value(s) - for s in start_value.free_symbols - {self.time} - if str(s) not in self.parameter_ids - } + changes += self._start( + period.time, self._with_model_values(start_value) ) - changes += self._start(period.time, start_value) elif ( {c.target_id for c in period_changes} & self.symbol_ids # preequilibration held the slope at zero or (self.ramps and periods[i_period - 1].is_preequilibration) + # a switch at the start of the period + or not (start_value - held).is_zero ): changes += self._start(period.time, start_value) for t, value in switch_values: changes += self._start(t, value) + held = previous return changes def _model_value(self, symbol: sp.Symbol) -> float: """The value of the model parameter ``symbol``.""" return self.sbml_model.getParameter(str(symbol)).getValue() + def _with_model_values(self, value: sp.Expr) -> sp.Expr: + """``value`` for a change at the start of an experiment, which may + only refer to parameter table symbols: the other symbols are constant + model parameters, at their values.""" + return value.subs( + { + s: self._model_value(s) + for s in value.free_symbols - {self.time} + if str(s) not in self.parameter_ids + } + ) + def _switching_times(self, subs: dict, experiment_id: str) -> list[float]: """The switching times, for the symbol values ``subs``.""" times = [] @@ -1109,6 +1139,41 @@ def _add_start_periods(experiments: list[Experiment]) -> set[int]: return droppable +def _referenced_ids( + problem: Problem, + additions: list[tuple[Experiment, float, str, sp.Expr]], +) -> set[str]: + """The IDs that the math of the model, the observables, or the condition + changes, including those in ``additions``, refer to.""" + ids = set() + for element in problem.model.sbml_model.getListOfAllElements(): + if not (hasattr(element, "isSetMath") and element.isSetMath()): + continue + nodes = [element.getMath()] + while nodes: + node = nodes.pop() + if node.getType() == libsbml.AST_NAME: + ids.add(node.getName()) + nodes += [node.getChild(i) for i in range(node.getNumChildren())] + expressions = [ + *(o.formula for o in problem.observables), + *(o.noise_formula for o in problem.observables), + *(value for *_, value in additions), + ] + ids |= { + str(s) + for expression in expressions + if expression is not None + for s in expression.free_symbols + } + ids |= { + str(s) + for table in problem.condition_tables + for s in table.free_symbols + } + return ids + + def _drop_changes( problem: Problem, experiments: list[Experiment], target_ids: set[str] ) -> None: diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index 201a6740..c558d9d7 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -670,7 +670,7 @@ def test_assignments_to_experiments(condition, before, after): species A, B; A = 10; B = 0; - k1 = piecewise(1, {condition}, 2); + k1 := piecewise(1, {condition}, 2); k2 = 0.5; R1: A -> B; k1 * A; R2: B -> A; k2 * B; @@ -719,7 +719,7 @@ def test_assignments_to_experiments_unused_condition_id(): problem = Problem() problem.model = SbmlModel.from_antimony(""" species A = 10; - k1 = piecewise(1, time < 10, 2); + k1 := piecewise(1, time < 10, 2); R1: A -> ; k1 * A; """) problem.add_observable("obs_a", formula="A", noise_formula="cond_0") @@ -1477,8 +1477,9 @@ def test_assignments_to_experiments_preequilibration_boundary( def test_assignments_to_experiments_validate(): """Test that the converted problem is linted. - Here, the input problem is valid, but the conversion makes the estimated - parameter `k1` a condition table target. + Here, the input problem is valid, but the conversion adds the model + parameter `k1_slope` for the ramp, which is also an estimated observable + parameter in the parameter table. """ problem = Problem() ant_model = """ @@ -1486,18 +1487,198 @@ def test_assignments_to_experiments_validate(): species A, B; A = 10; B = 0; - k1 = piecewise(0, time < 10, 1); + k1 := piecewise(time, time < 10, 1); k2 = 0.5; R1: A -> B; k1 * A; R2: B -> A; k2 * B; end """ problem.model = SbmlModel.from_antimony(ant_model) - problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_observable("obs_a", formula="k1_slope * A", noise_formula=1.0) problem.add_measurement("obs_a", time=20.0, measurement=5.0) - problem.add_parameter("k1", estimate=True, lb=1e-3, ub=1e3) + problem.add_parameter("k1_slope", estimate=True, lb=1e-3, ub=1e3) problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) assert not problem.validate().has_errors() with pytest.raises(ValueError, match="did not pass linting"): assignments_to_experiments(problem) + + +def _changes_by_time( + problem: Problem, experiment_id: str +) -> dict[float, dict[str, float]]: + """The changes of each period of an experiment, by the period's time.""" + return { + period.time: { + c.target_id: float(c.target_value) + for c in problem.get_changes_for_period(period) + } + for period in problem[experiment_id].sorted_periods + } + + +def _switch_model(assignment: str, extra: str = "") -> SbmlModel: + """A model with the assignment rule `k1 := `.""" + return SbmlModel.from_antimony(f""" + model conversion + species A; + A = 10; + other = 1; + {extra} + k1 := {assignment}; + R1: A -> ; k1 * A; + end + """) + + +def test_assignments_to_experiments_switch_at_period_start(): + """Test that a switch at the start of a later period is kept, although + that period does not change any symbol of the assignment.""" + problem = Problem() + problem.model = _switch_model("piecewise(0, time < 10, 1)") + problem.add_condition("c0", other=1.0) + problem.add_condition("c1", other=2.0) + problem.add_experiment("e1", 0, "c0", 10, "c1") + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + converted = assignments_to_experiments(problem) + + assert _changes_by_time(converted, "e1") == { + 0.0: {"other": 1.0}, + 10.0: {"other": 2.0, "k1": 1.0}, + } + + +def test_assignments_to_experiments_switch_at_preequilibration_end(): + """Test that a switch at time 0 is kept after preequilibration, where the + condition still holds at time 0.""" + problem = Problem() + problem.model = _switch_model("piecewise(0, time <= 0, 1)") + problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, [], 0, []) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + with pytest.warns(UserWarning, match="during preequilibration"): + converted = assignments_to_experiments(problem) + + assert _changes_by_time(converted, "e1") == { + v2.C.TIME_PREEQUILIBRATION: {"k1": 0.0}, + 0.0: {"k1": 1.0}, + } + + +def test_assignments_to_experiments_shared_symbol(): + """Test that the changes of a symbol of a converted assignment stay if + the model still refers to that symbol elsewhere.""" + problem = Problem() + problem.model = _switch_model( + "dose * piecewise(0, time < 10, 1)", + extra="dose = 0; R2: -> A; dose;", + ) + problem.add_condition("c1", dose=2.0) + problem.add_experiment("e1", 0, "c1") + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + converted = assignments_to_experiments(problem) + + assert _changes_by_time(converted, "e1") == { + 0.0: {"dose": 2.0}, + 10.0: {"k1": 2.0}, + } + + +def test_assignments_to_experiments_initial_assignment(): + """Test that a time-dependent initial assignment stays in the model, since + it is evaluated once, at the start of the simulation.""" + problem = Problem() + problem.model = SbmlModel.from_antimony(""" + model conversion + species A; + A = 10; + k1 = piecewise(0, time < 10, 1); + R1: A -> ; k1 * A; + end + """) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement("obs_a", time=20.0, measurement=1.0) + + converted = assignments_to_experiments(problem) + + assert converted.model.sbml_model.getInitialAssignment("k1") is not None + assert not converted.experiments + assert not converted.conditions + + +def test_assignments_to_experiments_preequilibration_model_value(): + """Test that a value at preequilibration, which starts the experiment, + only refers to parameter table symbols: other model parameters take their + values.""" + problem = Problem() + problem.model = _switch_model( + "piecewise(level, time < 10, 0)", extra="level = 3;" + ) + problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, [], 0, []) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + with pytest.warns(UserWarning, match="during preequilibration"): + converted = assignments_to_experiments(problem) + + assert _changes_by_time(converted, "e1") == { + v2.C.TIME_PREEQUILIBRATION: {"k1": 3.0}, + 0.0: {}, + 10.0: {"k1": 0.0}, + } + + +def test_assignments_to_experiments_species_target(): + """Test that an assignment to a species is rejected.""" + problem = Problem() + problem.model = SbmlModel.from_antimony(""" + model conversion + species A, S; + A = 10; + S := piecewise(0, time < 10, 1); + R1: A -> ; S * A; + end + """) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement("obs_a", time=20.0, measurement=1.0) + + with pytest.raises(NotImplementedError, match="not a parameter"): + assignments_to_experiments(problem) + + +def test_assignments_to_experiments_nothing_to_convert(): + """Test that a problem without time-dependent assignments is unchanged, + including an experiment that only preequilibrates.""" + problem = Problem() + problem.model = SbmlModel.from_antimony(""" + model conversion + species A; + A = 10; + k2 = 0.5; + R1: A -> ; k2 * A; + end + """) + problem.add_condition("c1", k2=1.0) + problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, "c1") + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + converted = assignments_to_experiments(problem) + + assert converted.experiment_df.equals(problem.experiment_df) + assert converted.condition_df.equals(problem.condition_df) From 520d565fcf93d364655eda7f65a884fe114c3740 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 16:10:25 +0000 Subject: [PATCH 13/15] assignments_to_experiments: count conversion factors as uses, leave unconverted problems alone, pick an unused slope ID - A parameter used as a model or species conversion factor counted as unused, so a condition change of it was dropped when the parameter also appeared in a converted assignment. Conversion factors now count as references. - A problem without time-dependent assignments was still changed: its experiments' periods were re-sorted. It is now returned as it was, after the lint. - A ramp's slope parameter was only named so as not to clash with SBML IDs, and clashed with, e.g., an observable parameter `_slope`. It is now named with `Problem.get_unused_id`. `test_assignments_to_experiments_validate` relied on that slope clash to fail the lint; it now uses an existing experiment `default_experiment`, which the measurements without an experiment ID are assigned to as well. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 33 +++++++++---- tests/v2/test_conversion.py | 96 ++++++++++++++++++++++++++++--------- 2 files changed, 98 insertions(+), 31 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index b0f5b0b5..dd14c366 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -820,6 +820,9 @@ def assignments_to_experiments(problem: Problem) -> Problem: expr = sbml_math_to_sympy(rule) if any(pw.has(TimeSymbol) for pw in expr.atoms(sp.Piecewise)): time_dependent.append((rule, expr)) + if not time_dependent: + _check_lint(problem) + return problem experiments = list(problem.experiments) if any(m.experiment_id is None for m in problem.measurements): @@ -827,7 +830,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: experiments.append(default_experiment) else: default_experiment = None - droppable = _add_start_periods(experiments) if time_dependent else set() + droppable = _add_start_periods(experiments) # (experiment, period start time, target ID, target value) additions = [] @@ -876,6 +879,13 @@ def assignments_to_experiments(problem: Problem) -> Problem: if measurement.experiment_id is None: measurement.experiment_id = DEFAULT_EXPERIMENT_ID + _check_lint(problem) + return problem + + +def _check_lint(problem: Problem) -> None: + """Raise a :class:`ValueError` if the converted ``problem`` does not + pass linting.""" validation_issues = problem.validate() if validation_issues.has_errors(): messages = "\n".join(f"{issue.message}" for issue in validation_issues) @@ -883,8 +893,6 @@ def assignments_to_experiments(problem: Problem) -> Problem: f"The converted PEtab problem did not pass linting:\n{messages}" ) - return problem - class _PiecewiseAssignment: """A time-dependent piecewise assignment to a model parameter, and the @@ -964,9 +972,9 @@ def replace_in_model(self) -> None: target's, is set by every period that starts one. """ if self.ramps: - self.slope_id = f"{self.target_id}_slope" - while self.sbml_model.getElementBySId(self.slope_id) is not None: - self.slope_id += "_" + self.slope_id = self.problem.get_unused_id( + f"{self.target_id}_slope" + ) add_sbml_parameter( self.sbml_model, self.slope_id, value=0.0, constant=False ) @@ -1143,10 +1151,11 @@ def _referenced_ids( problem: Problem, additions: list[tuple[Experiment, float, str, sp.Expr]], ) -> set[str]: - """The IDs that the math of the model, the observables, or the condition - changes, including those in ``additions``, refer to.""" + """The IDs that the model, the observables, or the condition changes, + including those in ``additions``, refer to.""" + sbml_model = problem.model.sbml_model ids = set() - for element in problem.model.sbml_model.getListOfAllElements(): + for element in sbml_model.getListOfAllElements(): if not (hasattr(element, "isSetMath") and element.isSetMath()): continue nodes = [element.getMath()] @@ -1155,6 +1164,12 @@ def _referenced_ids( if node.getType() == libsbml.AST_NAME: ids.add(node.getName()) nodes += [node.getChild(i) for i in range(node.getNumChildren())] + # conversion factors refer to parameters by attribute, not in math + ids |= { + element.getConversionFactor() + for element in [sbml_model, *sbml_model.getListOfSpecies()] + if element.isSetConversionFactor() + } expressions = [ *(o.formula for o in problem.observables), *(o.noise_formula for o in problem.observables), diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index c558d9d7..beacc493 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -1477,27 +1477,22 @@ def test_assignments_to_experiments_preequilibration_boundary( def test_assignments_to_experiments_validate(): """Test that the converted problem is linted. - Here, the input problem is valid, but the conversion adds the model - parameter `k1_slope` for the ramp, which is also an estimated observable - parameter in the parameter table. + Here, the input problem is valid, but the measurement without an + experiment ID is assigned to a new experiment `default_experiment`, an ID + that an existing experiment already has. """ problem = Problem() - ant_model = """ - model conversion - species A, B; - A = 10; - B = 0; - k1 := piecewise(time, time < 10, 1); - k2 = 0.5; - R1: A -> B; k1 * A; - R2: B -> A; k2 * B; - end - """ - problem.model = SbmlModel.from_antimony(ant_model) - problem.add_observable("obs_a", formula="k1_slope * A", noise_formula=1.0) - problem.add_measurement("obs_a", time=20.0, measurement=5.0) - problem.add_parameter("k1_slope", estimate=True, lb=1e-3, ub=1e3) - problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + problem.model = _switch_model("piecewise(0, time < 10, 1)") + problem.add_condition("c1", other=2.0) + problem.add_experiment(DEFAULT_EXPERIMENT_ID, 0, "c1") + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", + experiment_id=DEFAULT_EXPERIMENT_ID, + time=20.0, + measurement=5.0, + ) + problem.add_measurement("obs_a", time=30.0, measurement=5.0) assert not problem.validate().has_errors() with pytest.raises(ValueError, match="did not pass linting"): @@ -1661,7 +1656,8 @@ def test_assignments_to_experiments_species_target(): def test_assignments_to_experiments_nothing_to_convert(): """Test that a problem without time-dependent assignments is unchanged, - including an experiment that only preequilibrates.""" + including an experiment that only preequilibrates and periods that are + not in time order.""" problem = Problem() problem.model = SbmlModel.from_antimony(""" model conversion @@ -1672,7 +1668,43 @@ def test_assignments_to_experiments_nothing_to_convert(): end """) problem.add_condition("c1", k2=1.0) + problem.add_condition("c2", k2=2.0) problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, "c1") + # periods out of time order + problem.add_experiment("e2", 10, "c2", 0, "c1") + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + for experiment_id in ("e1", "e2"): + problem.add_measurement( + "obs_a", experiment_id=experiment_id, time=20.0, measurement=1.0 + ) + + converted = assignments_to_experiments(problem) + + def periods(problem): + return [ + [(p.time, p.condition_ids) for p in e.periods] + for e in problem.experiments + ] + + assert periods(converted) == periods(problem) + assert converted.condition_df.equals(problem.condition_df) + + +@pytest.mark.parametrize("on", ["model", "species"]) +def test_assignments_to_experiments_conversion_factor(on): + """Test that the changes of a symbol of a converted assignment stay if + the model refers to that symbol as a conversion factor.""" + problem = Problem() + problem.model = _switch_model( + "cf * piecewise(0, time < 10, 1)", extra="cf = 1;" + ) + sbml_model = problem.model.sbml_model + element = sbml_model if on == "model" else sbml_model.getSpecies("A") + assert ( + element.setConversionFactor("cf") == libsbml.LIBSBML_OPERATION_SUCCESS + ) + problem.add_condition("c1", cf=2.0) + problem.add_experiment("e1", 0, "c1") problem.add_observable("obs_a", formula="A", noise_formula=1.0) problem.add_measurement( "obs_a", experiment_id="e1", time=20.0, measurement=1.0 @@ -1680,5 +1712,25 @@ def test_assignments_to_experiments_nothing_to_convert(): converted = assignments_to_experiments(problem) - assert converted.experiment_df.equals(problem.experiment_df) - assert converted.condition_df.equals(problem.condition_df) + assert _changes_by_time(converted, "e1") == { + 0.0: {"cf": 2.0}, + 10.0: {"k1": 2.0}, + } + + +def test_assignments_to_experiments_slope_id_taken(): + """Test that a ramp's slope parameter gets an ID that the problem does not + use yet, here because an observable parameter is called `k1_slope`.""" + problem = Problem() + problem.model = _switch_model("piecewise(time, time < 10, 1)") + problem.add_observable("obs_a", formula="k1_slope * A", noise_formula=1.0) + problem.add_measurement("obs_a", time=20.0, measurement=5.0) + problem.add_parameter("k1_slope", estimate=True, lb=1e-3, ub=1e3) + + converted = assignments_to_experiments(problem) + + assert _rate_rule(converted, "k1") == "k1_slope_1" + assert _changes_by_time(converted, DEFAULT_EXPERIMENT_ID) == { + 0.0: {"k1": 0.0, "k1_slope_1": 1.0}, + 10.0: {"k1": 1.0, "k1_slope_1": 0.0}, + } From d9182ded882cccc9dff5a9065524b9cd554af450 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 16:21:52 +0000 Subject: [PATCH 14/15] assignments_to_experiments: evaluate a condition change once, at the start of its period The converter kept a condition change's value as an expression and substituted it wherever a switch used it, so the value was re-evaluated at the switch with the values at that time. A condition change is evaluated once, at the start of its period, and all changes of a period simultaneously, with the values from before the period. E.g., with `d1 = b1 + 0.5` at t = 2 and `b1 = 2` at t = 3, a switch to `d1` at t = 5 used 2.5 instead of 1.7; likewise with `b1 = 2` in the same condition as `d1`. A change's value is now evaluated at the start of its period: symbols set by an earlier condition take that value, other constant model parameters their model values, and parameter table symbols stay. A change of a symbol of the assignment whose value still refers to anything else, e.g. a species, is rejected, since its value is not known there. Found by differential testing against simulations of the original and the converted problem. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 45 ++++++++++++++++++++++++++++- tests/v2/test_conversion.py | 57 +++++++++++++++++++++++++++++++++++++ 2 files changed, 101 insertions(+), 1 deletion(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index dd14c366..b5666f04 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -1004,12 +1004,32 @@ def changes( held = None for i_period, period in enumerate(periods): period_changes = self.problem.get_changes_for_period(period) - overrides |= {c.target_id: c.target_value for c in period_changes} + # a condition's changes are evaluated once, at the start of their + # period, all with the values from before the period + overrides |= { + c.target_id: self._at_period_start(c.target_value, overrides) + for c in period_changes + } subs = self.fixed_values | { s: overrides[str(s)] for s in self.expr.free_symbols if str(s) in overrides } + if unknown := sorted( + { + str(x) + for s in self.expr.free_symbols + if str(s) in overrides + for x in overrides[str(s)].free_symbols + } + - self.parameter_ids + ): + raise NotImplementedError( + "Only condition changes whose values are known at the " + "start of their period are supported, but in experiment " + f"`{experiment.id}`, the changes of the symbols of the " + f"assignment to `{self.target_id}` refer to {unknown}." + ) times = self._switching_times(subs, experiment.id) if period.is_preequilibration: @@ -1075,6 +1095,29 @@ def _model_value(self, symbol: sp.Symbol) -> float: """The value of the model parameter ``symbol``.""" return self.sbml_model.getParameter(str(symbol)).getValue() + def _at_period_start(self, value: sp.Expr, overrides: dict) -> sp.Expr: + """The ``value`` of a change, evaluated at the start of its period. + + Symbols that a condition set before take that value (``overrides``), + and the other constant model parameters outside the parameter table + their model values. Parameter table symbols stay, since no condition + sets them; any other symbol stays as well, since its value is not + known here. + """ + if value is None: + return value + known = {} + for s in value.free_symbols: + if str(s) in overrides: + known[s] = overrides[str(s)] + elif ( + str(s) not in self.parameter_ids + and (p := self.sbml_model.getParameter(str(s))) is not None + and p.getConstant() + ): + known[s] = sp.Float(p.getValue()) + return value.xreplace(known) + def _with_model_values(self, value: sp.Expr) -> sp.Expr: """``value`` for a change at the start of an experiment, which may only refer to parameter table symbols: the other symbols are constant diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index beacc493..a3994cd4 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -1734,3 +1734,60 @@ def test_assignments_to_experiments_slope_id_taken(): 0.0: {"k1": 0.0, "k1_slope_1": 1.0}, 10.0: {"k1": 1.0, "k1_slope_1": 0.0}, } + + +@pytest.mark.parametrize( + "conditions", + [ + # the changes of a condition are simultaneous: `d1` takes the old `b1` + [(2.0, {"b1": 2.0, "d1": "b1 + 0.5"})], + # a change is evaluated once: `d1` keeps the `b1` of its period + [(2.0, {"d1": "b1 + 0.5"}), (3.0, {"b1": 2.0})], + ], +) +def test_assignments_to_experiments_change_evaluated_at_period_start( + conditions, +): + """Test that a condition change's value is evaluated at the start of its + period, with the values from before the period, also when a later switch + uses it.""" + problem = Problem() + problem.model = _switch_model( + "piecewise(b0, time < 5, d1)", extra="b0 = 0.3; b1 = 1.2; d1 = 1;" + ) + args = [0, []] + for i, (t, changes) in enumerate(conditions): + problem.add_condition(f"c{i}", **changes) + args += [t, f"c{i}"] + problem.add_experiment("e1", *args) + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + converted = assignments_to_experiments(problem) + + (switch,) = (p for p in converted["e1"].periods if p.time == 5.0) + (change,) = converted.get_changes_for_period(switch) + assert change.target_id == "k1" + # 1.2 + 0.5, not 2.0 + 0.5 + assert float(change.target_value) == pytest.approx(1.7) + + +def test_assignments_to_experiments_change_refers_to_species(): + """Test that a change of a symbol of the assignment whose value refers to + a species, which is not known at the start of the period, is + rejected.""" + problem = Problem() + problem.model = _switch_model( + "piecewise(0, time < 5, d1)", extra="d1 = 1;" + ) + problem.add_condition("c1", d1="A") + problem.add_experiment("e1", 0, [], 2, "c1") + problem.add_observable("obs_a", formula="A", noise_formula=1.0) + problem.add_measurement( + "obs_a", experiment_id="e1", time=20.0, measurement=1.0 + ) + + with pytest.raises(NotImplementedError, match=r"refer to \['A'\]"): + assignments_to_experiments(problem) From d9ac2965a165ed7021d14f4cfb9313688bc43740 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fabian=20Fr=C3=B6hlich?= Date: Sun, 25 Oct 2026 16:47:25 +0000 Subject: [PATCH 15/15] Move the fixes to their own PR This PR keeps only the refactor of `assignments_to_experiments`, which does not change its behaviour. The fixes, which do, follow in a separate PR stacked on this one. Co-Authored-By: Claude Opus 5.5 --- petab/v2/converters.py | 183 ++++----------------- tests/v2/test_conversion.py | 320 ++---------------------------------- 2 files changed, 45 insertions(+), 458 deletions(-) diff --git a/petab/v2/converters.py b/petab/v2/converters.py index b5666f04..9eef1061 100644 --- a/petab/v2/converters.py +++ b/petab/v2/converters.py @@ -757,9 +757,9 @@ def assignments_to_experiments(problem: Problem) -> Problem: if all of them were replaced, since a period's time overrides the model's initial time. - **Scope:** conversion is supported for SBML assignment rules to - parameters that depend on time through the conditions of piecewise - functions, and at most linearly through their values, where + **Scope:** conversion is supported for SBML assignment rules that depend + on time through the conditions of piecewise functions, and at most + linearly through their values, where * every such condition is an inequality that can be re-arranged explicitly into `t />= ...`, @@ -781,11 +781,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: entity instead, with the rate rule `d target/dt = _slope`, and each period sets both the target and its slope; an assignment that only switches between constant values stays a constant parameter. Anything else - raises a :class:`NotImplementedError`. Initial assignments stay in the - model, since they are evaluated once, at the start of the simulation. - - Condition changes of the replaced assignments' symbols are dropped, unless - the model, an observable or another change still refers to them. + raises a :class:`NotImplementedError`. .. warning:: @@ -811,19 +807,6 @@ def assignments_to_experiments(problem: Problem) -> Problem: problem = deepcopy(problem) sbml_model = problem.model.sbml_model - # Initial assignments are evaluated once, at the start of the - # simulation, so they never switch and stay in the model. - time_dependent = [] - for rule in sbml_model.getListOfRules(): - if not rule.isAssignment(): - continue - expr = sbml_math_to_sympy(rule) - if any(pw.has(TimeSymbol) for pw in expr.atoms(sp.Piecewise)): - time_dependent.append((rule, expr)) - if not time_dependent: - _check_lint(problem) - return problem - experiments = list(problem.experiments) if any(m.experiment_id is None for m in problem.measurements): default_experiment = Experiment(id=DEFAULT_EXPERIMENT_ID) @@ -836,7 +819,13 @@ def assignments_to_experiments(problem: Problem) -> Problem: additions = [] replaced_symbols = set() converted_ids = [] - for rule, expr in time_dependent: + for rule in [ + *(r for r in sbml_model.getListOfRules() if r.isAssignment()), + *sbml_model.getListOfInitialAssignments(), + ]: + expr = sbml_math_to_sympy(rule) + if not any(pw.has(TimeSymbol) for pw in expr.atoms(sp.Piecewise)): + continue assignment = _PiecewiseAssignment(rule, expr, problem) assignment.replace_in_model() for experiment in experiments: @@ -857,9 +846,7 @@ def assignments_to_experiments(problem: Problem) -> Problem: stacklevel=2, ) - # replace the changes that are covered by the new conditions, unless - # something else still refers to their targets - replaced_symbols -= _referenced_ids(problem, additions) + # replace the changes that are covered by the new conditions _drop_changes(problem, experiments, replaced_symbols) _add_conditions(problem, additions) @@ -879,13 +866,6 @@ def assignments_to_experiments(problem: Problem) -> Problem: if measurement.experiment_id is None: measurement.experiment_id = DEFAULT_EXPERIMENT_ID - _check_lint(problem) - return problem - - -def _check_lint(problem: Problem) -> None: - """Raise a :class:`ValueError` if the converted ``problem`` does not - pass linting.""" validation_issues = problem.validate() if validation_issues.has_errors(): messages = "\n".join(f"{issue.message}" for issue in validation_issues) @@ -893,6 +873,8 @@ def _check_lint(problem: Problem) -> None: f"The converted PEtab problem did not pass linting:\n{messages}" ) + return problem + class _PiecewiseAssignment: """A time-dependent piecewise assignment to a model parameter, and the @@ -911,11 +893,6 @@ def __init__( self.problem = problem self.sbml_model = problem.model.sbml_model self.target_id = rule.getId() - if self.sbml_model.getParameter(self.target_id) is None: - raise NotImplementedError( - "Only assignments to parameters are supported, but " - f"`{self.target_id}` is not a parameter." - ) self.expr = expr (self.time,) = expr.atoms(TimeSymbol) self.symbol_ids = {str(s) for s in expr.free_symbols} @@ -972,9 +949,9 @@ def replace_in_model(self) -> None: target's, is set by every period that starts one. """ if self.ramps: - self.slope_id = self.problem.get_unused_id( - f"{self.target_id}_slope" - ) + self.slope_id = f"{self.target_id}_slope" + while self.sbml_model.getElementBySId(self.slope_id) is not None: + self.slope_id += "_" add_sbml_parameter( self.sbml_model, self.slope_id, value=0.0, constant=False ) @@ -1000,50 +977,22 @@ def changes( changes = [] periods = experiment.sorted_periods overrides = {} - # the value at the end of the previous period - held = None for i_period, period in enumerate(periods): period_changes = self.problem.get_changes_for_period(period) - # a condition's changes are evaluated once, at the start of their - # period, all with the values from before the period - overrides |= { - c.target_id: self._at_period_start(c.target_value, overrides) - for c in period_changes - } + overrides |= {c.target_id: c.target_value for c in period_changes} subs = self.fixed_values | { s: overrides[str(s)] for s in self.expr.free_symbols if str(s) in overrides } - if unknown := sorted( - { - str(x) - for s in self.expr.free_symbols - if str(s) in overrides - for x in overrides[str(s)].free_symbols - } - - self.parameter_ids - ): - raise NotImplementedError( - "Only condition changes whose values are known at the " - "start of their period are supported, but in experiment " - f"`{experiment.id}`, the changes of the symbols of the " - f"assignment to `{self.target_id}` refer to {unknown}." - ) times = self._switching_times(subs, experiment.id) if period.is_preequilibration: # the value at time 0, which respects whether a condition # holds for a switch at time 0, held there: a ramp does not # run during preequilibration - held = self.expr.subs(subs).subs(self.time, 0) - changes.append( - ( - period.time, - self.target_id, - self._with_model_values(held), - ) - ) + value = self.expr.subs(subs).subs(self.time, 0) + changes.append((period.time, self.target_id, value)) if self.ramps: changes.append((period.time, self.slope_id, sp.S.Zero)) continue @@ -1075,61 +1024,31 @@ def changes( or {str(s) for s in self.default.free_symbols} & self.parameter_ids ): - changes += self._start( - period.time, self._with_model_values(start_value) + # at the start of an experiment, a change may only refer + # to parameter table symbols; other symbols are + # constant model parameters, at their values + start_value = start_value.subs( + { + s: self._model_value(s) + for s in start_value.free_symbols - {self.time} + if str(s) not in self.parameter_ids + } ) + changes += self._start(period.time, start_value) elif ( {c.target_id for c in period_changes} & self.symbol_ids # preequilibration held the slope at zero or (self.ramps and periods[i_period - 1].is_preequilibration) - # a switch at the start of the period - or not (start_value - held).is_zero ): changes += self._start(period.time, start_value) for t, value in switch_values: changes += self._start(t, value) - held = previous return changes def _model_value(self, symbol: sp.Symbol) -> float: """The value of the model parameter ``symbol``.""" return self.sbml_model.getParameter(str(symbol)).getValue() - def _at_period_start(self, value: sp.Expr, overrides: dict) -> sp.Expr: - """The ``value`` of a change, evaluated at the start of its period. - - Symbols that a condition set before take that value (``overrides``), - and the other constant model parameters outside the parameter table - their model values. Parameter table symbols stay, since no condition - sets them; any other symbol stays as well, since its value is not - known here. - """ - if value is None: - return value - known = {} - for s in value.free_symbols: - if str(s) in overrides: - known[s] = overrides[str(s)] - elif ( - str(s) not in self.parameter_ids - and (p := self.sbml_model.getParameter(str(s))) is not None - and p.getConstant() - ): - known[s] = sp.Float(p.getValue()) - return value.xreplace(known) - - def _with_model_values(self, value: sp.Expr) -> sp.Expr: - """``value`` for a change at the start of an experiment, which may - only refer to parameter table symbols: the other symbols are constant - model parameters, at their values.""" - return value.subs( - { - s: self._model_value(s) - for s in value.free_symbols - {self.time} - if str(s) not in self.parameter_ids - } - ) - def _switching_times(self, subs: dict, experiment_id: str) -> list[float]: """The switching times, for the symbol values ``subs``.""" times = [] @@ -1190,48 +1109,6 @@ def _add_start_periods(experiments: list[Experiment]) -> set[int]: return droppable -def _referenced_ids( - problem: Problem, - additions: list[tuple[Experiment, float, str, sp.Expr]], -) -> set[str]: - """The IDs that the model, the observables, or the condition changes, - including those in ``additions``, refer to.""" - sbml_model = problem.model.sbml_model - ids = set() - for element in sbml_model.getListOfAllElements(): - if not (hasattr(element, "isSetMath") and element.isSetMath()): - continue - nodes = [element.getMath()] - while nodes: - node = nodes.pop() - if node.getType() == libsbml.AST_NAME: - ids.add(node.getName()) - nodes += [node.getChild(i) for i in range(node.getNumChildren())] - # conversion factors refer to parameters by attribute, not in math - ids |= { - element.getConversionFactor() - for element in [sbml_model, *sbml_model.getListOfSpecies()] - if element.isSetConversionFactor() - } - expressions = [ - *(o.formula for o in problem.observables), - *(o.noise_formula for o in problem.observables), - *(value for *_, value in additions), - ] - ids |= { - str(s) - for expression in expressions - if expression is not None - for s in expression.free_symbols - } - ids |= { - str(s) - for table in problem.condition_tables - for s in table.free_symbols - } - return ids - - def _drop_changes( problem: Problem, experiments: list[Experiment], target_ids: set[str] ) -> None: diff --git a/tests/v2/test_conversion.py b/tests/v2/test_conversion.py index a3994cd4..201a6740 100644 --- a/tests/v2/test_conversion.py +++ b/tests/v2/test_conversion.py @@ -670,7 +670,7 @@ def test_assignments_to_experiments(condition, before, after): species A, B; A = 10; B = 0; - k1 := piecewise(1, {condition}, 2); + k1 = piecewise(1, {condition}, 2); k2 = 0.5; R1: A -> B; k1 * A; R2: B -> A; k2 * B; @@ -719,7 +719,7 @@ def test_assignments_to_experiments_unused_condition_id(): problem = Problem() problem.model = SbmlModel.from_antimony(""" species A = 10; - k1 := piecewise(1, time < 10, 2); + k1 = piecewise(1, time < 10, 2); R1: A -> ; k1 * A; """) problem.add_observable("obs_a", formula="A", noise_formula="cond_0") @@ -1477,317 +1477,27 @@ def test_assignments_to_experiments_preequilibration_boundary( def test_assignments_to_experiments_validate(): """Test that the converted problem is linted. - Here, the input problem is valid, but the measurement without an - experiment ID is assigned to a new experiment `default_experiment`, an ID - that an existing experiment already has. + Here, the input problem is valid, but the conversion makes the estimated + parameter `k1` a condition table target. """ problem = Problem() - problem.model = _switch_model("piecewise(0, time < 10, 1)") - problem.add_condition("c1", other=2.0) - problem.add_experiment(DEFAULT_EXPERIMENT_ID, 0, "c1") - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", - experiment_id=DEFAULT_EXPERIMENT_ID, - time=20.0, - measurement=5.0, - ) - problem.add_measurement("obs_a", time=30.0, measurement=5.0) - assert not problem.validate().has_errors() - - with pytest.raises(ValueError, match="did not pass linting"): - assignments_to_experiments(problem) - - -def _changes_by_time( - problem: Problem, experiment_id: str -) -> dict[float, dict[str, float]]: - """The changes of each period of an experiment, by the period's time.""" - return { - period.time: { - c.target_id: float(c.target_value) - for c in problem.get_changes_for_period(period) - } - for period in problem[experiment_id].sorted_periods - } - - -def _switch_model(assignment: str, extra: str = "") -> SbmlModel: - """A model with the assignment rule `k1 := `.""" - return SbmlModel.from_antimony(f""" - model conversion - species A; - A = 10; - other = 1; - {extra} - k1 := {assignment}; - R1: A -> ; k1 * A; - end - """) - - -def test_assignments_to_experiments_switch_at_period_start(): - """Test that a switch at the start of a later period is kept, although - that period does not change any symbol of the assignment.""" - problem = Problem() - problem.model = _switch_model("piecewise(0, time < 10, 1)") - problem.add_condition("c0", other=1.0) - problem.add_condition("c1", other=2.0) - problem.add_experiment("e1", 0, "c0", 10, "c1") - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) - - converted = assignments_to_experiments(problem) - - assert _changes_by_time(converted, "e1") == { - 0.0: {"other": 1.0}, - 10.0: {"other": 2.0, "k1": 1.0}, - } - - -def test_assignments_to_experiments_switch_at_preequilibration_end(): - """Test that a switch at time 0 is kept after preequilibration, where the - condition still holds at time 0.""" - problem = Problem() - problem.model = _switch_model("piecewise(0, time <= 0, 1)") - problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, [], 0, []) - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) - - with pytest.warns(UserWarning, match="during preequilibration"): - converted = assignments_to_experiments(problem) - - assert _changes_by_time(converted, "e1") == { - v2.C.TIME_PREEQUILIBRATION: {"k1": 0.0}, - 0.0: {"k1": 1.0}, - } - - -def test_assignments_to_experiments_shared_symbol(): - """Test that the changes of a symbol of a converted assignment stay if - the model still refers to that symbol elsewhere.""" - problem = Problem() - problem.model = _switch_model( - "dose * piecewise(0, time < 10, 1)", - extra="dose = 0; R2: -> A; dose;", - ) - problem.add_condition("c1", dose=2.0) - problem.add_experiment("e1", 0, "c1") - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) - - converted = assignments_to_experiments(problem) - - assert _changes_by_time(converted, "e1") == { - 0.0: {"dose": 2.0}, - 10.0: {"k1": 2.0}, - } - - -def test_assignments_to_experiments_initial_assignment(): - """Test that a time-dependent initial assignment stays in the model, since - it is evaluated once, at the start of the simulation.""" - problem = Problem() - problem.model = SbmlModel.from_antimony(""" + ant_model = """ model conversion - species A; + species A, B; A = 10; + B = 0; k1 = piecewise(0, time < 10, 1); - R1: A -> ; k1 * A; - end - """) - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement("obs_a", time=20.0, measurement=1.0) - - converted = assignments_to_experiments(problem) - - assert converted.model.sbml_model.getInitialAssignment("k1") is not None - assert not converted.experiments - assert not converted.conditions - - -def test_assignments_to_experiments_preequilibration_model_value(): - """Test that a value at preequilibration, which starts the experiment, - only refers to parameter table symbols: other model parameters take their - values.""" - problem = Problem() - problem.model = _switch_model( - "piecewise(level, time < 10, 0)", extra="level = 3;" - ) - problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, [], 0, []) - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) - - with pytest.warns(UserWarning, match="during preequilibration"): - converted = assignments_to_experiments(problem) - - assert _changes_by_time(converted, "e1") == { - v2.C.TIME_PREEQUILIBRATION: {"k1": 3.0}, - 0.0: {}, - 10.0: {"k1": 0.0}, - } - - -def test_assignments_to_experiments_species_target(): - """Test that an assignment to a species is rejected.""" - problem = Problem() - problem.model = SbmlModel.from_antimony(""" - model conversion - species A, S; - A = 10; - S := piecewise(0, time < 10, 1); - R1: A -> ; S * A; - end - """) - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement("obs_a", time=20.0, measurement=1.0) - - with pytest.raises(NotImplementedError, match="not a parameter"): - assignments_to_experiments(problem) - - -def test_assignments_to_experiments_nothing_to_convert(): - """Test that a problem without time-dependent assignments is unchanged, - including an experiment that only preequilibrates and periods that are - not in time order.""" - problem = Problem() - problem.model = SbmlModel.from_antimony(""" - model conversion - species A; - A = 10; k2 = 0.5; - R1: A -> ; k2 * A; + R1: A -> B; k1 * A; + R2: B -> A; k2 * B; end - """) - problem.add_condition("c1", k2=1.0) - problem.add_condition("c2", k2=2.0) - problem.add_experiment("e1", v2.C.TIME_PREEQUILIBRATION, "c1") - # periods out of time order - problem.add_experiment("e2", 10, "c2", 0, "c1") - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - for experiment_id in ("e1", "e2"): - problem.add_measurement( - "obs_a", experiment_id=experiment_id, time=20.0, measurement=1.0 - ) - - converted = assignments_to_experiments(problem) - - def periods(problem): - return [ - [(p.time, p.condition_ids) for p in e.periods] - for e in problem.experiments - ] - - assert periods(converted) == periods(problem) - assert converted.condition_df.equals(problem.condition_df) - - -@pytest.mark.parametrize("on", ["model", "species"]) -def test_assignments_to_experiments_conversion_factor(on): - """Test that the changes of a symbol of a converted assignment stay if - the model refers to that symbol as a conversion factor.""" - problem = Problem() - problem.model = _switch_model( - "cf * piecewise(0, time < 10, 1)", extra="cf = 1;" - ) - sbml_model = problem.model.sbml_model - element = sbml_model if on == "model" else sbml_model.getSpecies("A") - assert ( - element.setConversionFactor("cf") == libsbml.LIBSBML_OPERATION_SUCCESS - ) - problem.add_condition("c1", cf=2.0) - problem.add_experiment("e1", 0, "c1") + """ + problem.model = SbmlModel.from_antimony(ant_model) problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) - - converted = assignments_to_experiments(problem) - - assert _changes_by_time(converted, "e1") == { - 0.0: {"cf": 2.0}, - 10.0: {"k1": 2.0}, - } - - -def test_assignments_to_experiments_slope_id_taken(): - """Test that a ramp's slope parameter gets an ID that the problem does not - use yet, here because an observable parameter is called `k1_slope`.""" - problem = Problem() - problem.model = _switch_model("piecewise(time, time < 10, 1)") - problem.add_observable("obs_a", formula="k1_slope * A", noise_formula=1.0) problem.add_measurement("obs_a", time=20.0, measurement=5.0) - problem.add_parameter("k1_slope", estimate=True, lb=1e-3, ub=1e3) - - converted = assignments_to_experiments(problem) - - assert _rate_rule(converted, "k1") == "k1_slope_1" - assert _changes_by_time(converted, DEFAULT_EXPERIMENT_ID) == { - 0.0: {"k1": 0.0, "k1_slope_1": 1.0}, - 10.0: {"k1": 1.0, "k1_slope_1": 0.0}, - } - - -@pytest.mark.parametrize( - "conditions", - [ - # the changes of a condition are simultaneous: `d1` takes the old `b1` - [(2.0, {"b1": 2.0, "d1": "b1 + 0.5"})], - # a change is evaluated once: `d1` keeps the `b1` of its period - [(2.0, {"d1": "b1 + 0.5"}), (3.0, {"b1": 2.0})], - ], -) -def test_assignments_to_experiments_change_evaluated_at_period_start( - conditions, -): - """Test that a condition change's value is evaluated at the start of its - period, with the values from before the period, also when a later switch - uses it.""" - problem = Problem() - problem.model = _switch_model( - "piecewise(b0, time < 5, d1)", extra="b0 = 0.3; b1 = 1.2; d1 = 1;" - ) - args = [0, []] - for i, (t, changes) in enumerate(conditions): - problem.add_condition(f"c{i}", **changes) - args += [t, f"c{i}"] - problem.add_experiment("e1", *args) - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) - - converted = assignments_to_experiments(problem) - - (switch,) = (p for p in converted["e1"].periods if p.time == 5.0) - (change,) = converted.get_changes_for_period(switch) - assert change.target_id == "k1" - # 1.2 + 0.5, not 2.0 + 0.5 - assert float(change.target_value) == pytest.approx(1.7) - - -def test_assignments_to_experiments_change_refers_to_species(): - """Test that a change of a symbol of the assignment whose value refers to - a species, which is not known at the start of the period, is - rejected.""" - problem = Problem() - problem.model = _switch_model( - "piecewise(0, time < 5, d1)", extra="d1 = 1;" - ) - problem.add_condition("c1", d1="A") - problem.add_experiment("e1", 0, [], 2, "c1") - problem.add_observable("obs_a", formula="A", noise_formula=1.0) - problem.add_measurement( - "obs_a", experiment_id="e1", time=20.0, measurement=1.0 - ) + problem.add_parameter("k1", estimate=True, lb=1e-3, ub=1e3) + problem.add_parameter("k2", estimate=True, lb=1e-3, ub=1e3) + assert not problem.validate().has_errors() - with pytest.raises(NotImplementedError, match=r"refer to \['A'\]"): + with pytest.raises(ValueError, match="did not pass linting"): assignments_to_experiments(problem)