Usage
ode2sbml converts Python ODE systems into SBML Level 3 Version 2 and PEtab bundles.
Use the Python API when you want explicit control over model construction, conversion, validation, and PEtab export. Use the CLI when you want quick file-based conversion with auto-detection.
Input conventions
ode2sbml supports three model input styles.
Convention A: scipy-style function
Use this when your ODE model is written as a Python function with signature f(t, y, p) or f(t, y).
from ode2sbml import from_function, to_sbml
def sir(t, y, p):
S, I, R = y
beta, gamma = p
dS = -beta * S * I
dI = beta * S * I - gamma * I
dR = gamma * I
return [dS, dI, dR]
ir = from_function(
func=sir,
state_vars=["S", "I", "R"],
param_names=["beta", "gamma"],
initial_conditions={
"S": 990.0,
"I": 10.0,
"R": 0.0,
},
parameter_values={
"beta": 0.3,
"gamma": 0.1,
},
model_id="sir_model",
model_name="SIR Epidemic Model",
notes="Classic SIR compartmental model for infectious disease.",
)
to_sbml(ir, "sir_model.xml")
Best for existing simulation code
Convention A is useful when you already have scipy-style ODE code and want to convert it without rewriting the model as a dictionary or symbolic system.
Convention B: dict-of-formulas
Use this when you want a compact, readable model specification.
from ode2sbml import from_dict, to_sbml
model = {
"species": {
"S": 10.0,
"P": 0.0,
},
"parameters": {
"k": 0.1,
},
"odes": {
"S": "-k*S",
"P": "k*S",
},
}
ir = from_dict(
model,
model_id="conversion",
model_name="S to P conversion",
)
to_sbml(ir, "conversion.xml")
A larger dict model can also include compartments, assignment rules, and events:
from ode2sbml import from_dict, to_sbml
spec = {
"species": {
"x": 10.0,
"y": 5.0,
},
"parameters": {
"alpha": 1.1,
"beta": 0.4,
"delta": 0.1,
"gamma": 0.4,
},
"odes": {
"x": "alpha*x - beta*x*y",
"y": "-gamma*y + delta*x*y",
},
"compartments": {
"env": 1.0,
},
"assignment_rules": {
"total": "x + y",
},
"events": [
{
"trigger": "t >= 10",
"assignments": {
"x": "x * 2",
},
"delay": 0.0,
}
],
}
ir = from_dict(
spec,
model_id="lv",
model_name="Lotka-Volterra",
)
to_sbml(ir, "lotka_volterra.xml")
Recommended starting point
Convention B is usually the easiest format for writing reproducible examples, tests, tutorials, and small mechanistic models.
Convention C: SymPy symbolic system
Use this when your model is already symbolic or when you want exact symbolic expressions before export.
import sympy as sp
from ode2sbml import from_sympy, to_sbml
x, y, z = sp.symbols("x y z")
sigma, rho, b = sp.symbols("sigma rho b")
ir = from_sympy(
odes={
x: sigma * (y - x),
y: x * (rho - z) - y,
z: x * y - b * z,
},
params={
sigma: 10.0,
rho: 28.0,
b: 8.0 / 3.0,
},
inits={
x: 0.0,
y: 1.0,
z: 0.0,
},
model_id="lorenz",
model_name="Lorenz Attractor",
)
to_sbml(ir, "lorenz.xml")
Use SymPy for algebraic models
Convention C is useful when the model is generated from symbolic manipulation, notebooks, or analytical derivations.
Validate generated SBML
Validate an SBML file with the validator module.
from ode2sbml.validators.sbml_validator import validate_sbml_file
messages = validate_sbml_file(
"conversion.xml",
raise_on_error=False,
)
for message in messages:
print(message)
To fail immediately on SBML errors:
from ode2sbml.validators.sbml_validator import validate_sbml_file
validate_sbml_file(
"conversion.xml",
raise_on_error=True,
)
Validation is not optional
Always validate generated SBML before using it in downstream simulators, PEtab workflows, or publication pipelines.
Convert to PEtab
Export a PEtab bundle from an ODEModel.
from ode2sbml import from_dict, to_petab
model = {
"species": {
"S": 10.0,
"P": 0.0,
},
"parameters": {
"k": 0.1,
},
"odes": {
"S": "-k*S",
"P": "k*S",
},
}
ir = from_dict(
model,
model_id="conversion",
model_name="S to P conversion",
)
yaml_path = to_petab(
ir,
"build/conversion_petab",
)
print(yaml_path)
With custom observables:
from ode2sbml import from_dict, to_petab
ir = from_dict(
spec,
model_id="lv",
model_name="Lotka-Volterra",
)
yaml_path = to_petab(
ir,
outdir="petab_bundle",
observables={
"obs_x": {
"formula": "x",
"noise_formula": "sigma_x",
},
"obs_y": {
"formula": "y",
"noise_formula": "sigma_y",
},
},
)
print(yaml_path)
The PEtab export writes a bundle containing:
petab_bundle/
├── model.xml
├── conditions.tsv
├── observables.tsv
├── measurements.tsv
├── parameters.tsv
└── petab_problem.yaml
Default observables
If no observables are supplied, ode2sbml creates one observable per species.
Work directly with the IR
The central intermediate representation is ODEModel.
Use it directly when you want maximum control or when testing exporters without parser-specific assumptions.
from ode2sbml.model import (
ODEModel,
Compartment,
Species,
Parameter,
RateRule,
)
ir = ODEModel(
id="decay",
name="Decay model",
compartments=[
Compartment(
id="cell",
size=1.0,
)
],
species=[
Species(
id="x",
compartment="cell",
initial_value=1.0,
)
],
parameters=[
Parameter(
id="k",
value=0.2,
)
],
rate_rules=[
RateRule(
variable="x",
formula="-k*x",
)
],
)
json_text = ir.to_json()
ir2 = ODEModel.from_json(json_text)
Use the IR for tests
Constructing ODEModel directly is the cleanest way to unit-test exporters because it bypasses parser-specific assumptions.
JSON serialisation
ODEModel is serialisable through Pydantic.
from ode2sbml.model import ODEModel
json_text = ir.to_json()
restored_ir = ODEModel.from_json(json_text)
This is useful for caching, debugging, regression tests, and separating parsing from export.
Formula rules
Formulas should be SymPy-parseable strings.
Valid examples:
Avoid Python-only prefixes:
SBML identifier constraint
Model, species, compartment, and parameter IDs must match ^[A-Za-z_][A-Za-z0-9_]*$.
CLI usage
The CLI can auto-detect the model convention from a Python file.
Export a PEtab bundle:
Validate SBML:
The CLI detects input conventions in this order:
- Module-level
modeldictionary with anodeskey → Convention B. - Module-level
odesdictionary with SymPy symbol keys → Convention C. - Callable named
model,ode,f, orrhs→ Convention A.
For Convention A through the CLI, the Python module must also define state_vars.
Example CLI-compatible file:
def model(t, y, p):
S, P = y
k = p[0]
dS = -k * S
dP = k * S
return [dS, dP]
state_vars = ["S", "P"]
param_names = ["k"]
initial_conditions = {
"S": 10.0,
"P": 0.0,
}
parameter_values = {
"k": 0.1,
}
Then run:
CLI for fast conversion
Use the CLI for quick model conversion. Use the Python API when you need custom observables, custom PEtab metadata, direct IR manipulation, or integration into a larger workflow.
Common workflow
A typical programmatic workflow looks like this:
from ode2sbml import from_dict, to_sbml, to_petab
from ode2sbml.validators.sbml_validator import validate_sbml_file
spec = {
"species": {
"S": 10.0,
"P": 0.0,
},
"parameters": {
"k": 0.1,
},
"odes": {
"S": "-k*S",
"P": "k*S",
},
}
ir = from_dict(
spec,
model_id="conversion",
model_name="S to P conversion",
)
to_sbml(ir, "conversion.xml")
messages = validate_sbml_file(
"conversion.xml",
raise_on_error=False,
)
for message in messages:
print(message)
yaml_path = to_petab(
ir,
"build/conversion_petab",
)
print(yaml_path)
Production recommendation
For reproducible workflows, keep your model specification, generated SBML, validation logs, and PEtab bundle under version control or workflow-managed output directories.