Skip to content

assignments_to_experiments: keep the period that starts an experiment #559

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We鈥檒l occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
FFroehlich wants to merge 1 commit into from
Draft

assignments_to_experiments: keep the period that starts an experiment #559

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
39 changes: 23 additions & 16 deletions petab/v2/converters.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -671,7 +671,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
Expand Down Expand Up @@ -735,13 +738,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 = []
Expand Down Expand Up @@ -872,21 +884,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
Expand Down Expand Up @@ -973,7 +980,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
Expand Down
128 changes: 104 additions & 24 deletions tests/v2/test_conversion.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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"],
}
)

Expand Down Expand Up @@ -739,21 +741,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"],
}
)

Expand Down Expand Up @@ -847,14 +850,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"],
}
)

Expand All @@ -875,6 +880,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."""
Expand Down Expand Up @@ -1022,13 +1097,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)]


Expand Down Expand Up @@ -1118,14 +1196,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


Expand Down
Loading