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 = <target>_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 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01MNJePorYs9zmPYMam6g4EP

Author: FFroehlich

petab/v2/converters.py

@@ -674,7 +674,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 </<=/>/>= ...`,
@@ -688,7 +689,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 = <target>_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::
@@ -767,10 +774,27 @@ 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}
         }
         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 = <target>_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}
@@ -802,11 +826,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
@@ -832,14 +861,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(
@@ -913,6 +955,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,
@@ -972,7 +1024,12 @@ 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}`."
+        )
     return switches

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
@@ -888,7 +889,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)", "evaluate to a real number"),
         ("piecewise(0, time < A, 1)", "constant model parameters"),
         ("A * piecewise(0, time < 10, 1)", "constant model parameters"),
@@ -1036,6 +1038,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",