Model Definition

The Model class serves as the foundational component for describing biochemical systems in Catalax. It provides a unified framework for defining states, such as molecules, proteins or process variables, as well as reactions, parameters, and constraints that govern system dynamics. Whether you’re modeling simple enzyme kinetics or complex metabolic networks, the Model class offers the flexibility and mathematical rigor needed for advanced computational analysis including simulation, optimization, and Bayesian inference.

Understanding the Model Architecture

Core Components and Structure

The Model class organizes biochemical systems into five main components, each serving a specific purpose in mathematical modeling:

• Species: Chemical entities whose concentrations change over time (substrates, products, enzymes, inhibitors)
• ODEs: Differential equations governing the rate of change for each state
• Parameters: Kinetic constants and other model coefficients that characterize system behavior
• Constants: Time-invariant quantities that influence dynamics but don’t change during simulation
• Assignments: Derived quantities calculated from other model components

This modular architecture enables systematic model construction while maintaining mathematical consistency and biological interpretability.

Mathematical Foundation

The Model class implements systems of ordinary differential equations (ODEs) of the form:

$$\frac{dy_i}{dt} = f_i(y, \theta, c, t)$$

where $y_i$ represents the state, $\theta$ contains parameters, $c$ represents constants, and $t$ is time. This mathematical framework supports a wide range of biochemical phenomena from simple mass action kinetics to complex regulatory mechanisms.

Basic Model Construction Workflow

Step 1: Model Instantiation

Begin by creating a model instance with a descriptive name that clearly identifies your biochemical system:

import catalax as ctx

# Create a model with descriptive naming
model = ctx.Model(name="Enzyme Kinetics Study")

The model name serves important documentation purposes and appears in output files, visualizations, and analysis reports. Choose names that will be meaningful in research publications and collaborative work.

Step 2: Species Definition

Add states that participate in your system using flexible syntax options:

# Method 1: Keyword arguments with descriptive names
model.add_states(
    S="Substrate",
    E="Enzyme", 
    P="Product",
    ES="Enzyme-Substrate Complex"
)

# Method 2: Dictionary unpacking for programmatic construction
state_dict = {"I": "Inhibitor", "EI": "Enzyme-Inhibitor Complex"}
model.add_states(**state_dict)

States naming considerations:

• Symbols: Short identifiers used in equations (S, E, P)
• Names: Descriptive labels for documentation and visualization
• Consistency: Maintain consistent naming conventions across related models
• Biological meaning: Choose symbols that reflect chemical reality and facilitate interpretation

Step 3a: Reaction Definition

Define the reactions that govern states dynamics by either by providing a reaction scheme or by defining the reactants, products individually:

# Method 1: By providing a reaction scheme
model.add_reaction(
    "S -> P",
    symbol="r1",
    equation="k_cat * E * S / (K_m + S)",
)

# Method 2: By defining the reactants and products individually
model.add_reaction(
    "S -> P",
    symbol="r1",
    reactants=[("E", 1), ("S", 1)],
    products=[("P", 1)],
)

Key features of reaction definition:

• Automatic parameter inference: Parameters (k_cat, K_m) are automatically detected and added to the model
• Symbolic processing: Equations are converted to SymPy expressions for mathematical manipulation
• Observable states: By default, all states are observable; use observable=False for hidden states. In the case of using the scheme method, simply change the arrow to a double arrow to imply reversibility. For instance, S -> P becomes S <-> P.
• Validation: The system checks for invalid states not yet added to the model.

Reaction scheme supported arrows:

• -> (irreversible)
• <-> (reversible)
• => (irreversible)
• <=> (reversible)
• <==> (reversible)
• <===> (reversible)
• > (irreversible)

Step 3b: ODE Definition

Define the differential equations that govern states dynamics by either using the add_ode method or the add_odes method:

# Method 1: Individual ODEs
model.add_ode("S", "-(k_cat * E * S) / (K_m + S)")
model.add_ode("E", "0")  # Enzyme conservation
model.add_ode("P", "(k_cat * E * S) / (K_m + S)")
model.add_ode("ES", "0")  # Assuming rapid equilibrium

# Method 2: Multiple ODEs at once
model.add_odes(
    S="-(k_cat * E * S) / (K_m + S)",
    E="0",
    P="(k_cat * E * S) / (K_m + S)",
    ES="0"
)

Key features of ODE definition:

• Automatic parameter inference: Parameters (k_cat, K_m) are automatically detected and added to the model
• Symbolic processing: Equations are converted to SymPy expressions for mathematical manipulation
• Observable states: By default, all states are observable; use observable=False for hidden states
• Validation: The system checks for invalid states not yet added to the model.

Step 4: Parameter Configuration

Configure model parameters with values, bounds, and constraints:

# Set parameter values
model.parameters["k_cat"].value = 7.0
model.parameters["K_m"].value = 100.0

# Set parameter bounds for optimization
model.parameters["k_cat"].lower_bound = 0.1
model.parameters["k_cat"].upper_bound = 50.0

# Add parameter descriptions for documentation
model.parameters["k_cat"].name = "Turnover number"
model.parameters["K_m"].name = "Michaelis constant"

Advanced Model Features

Constants and External Parameters

Define time-invariant quantities that influence system behavior by either using the add_constant method or the add_constants method:

# Add constants that remain fixed during simulation
model.add_constant("pH", "solution_pH")
model.add_constant("temperature", "reaction_temperature") 

# Add constants that remain fixed during simulation
model.add_constants(
    pH="solution_pH",
    temperature="reaction_temperature"
)

# Constants can be accessed and used in ODEs
model.add_ode("S", "-(k_cat * E * S * pH) / (K_m + S)")

Constants are particularly useful for:

• Environmental conditions: pH, temperature, ionic strength
• Total concentrations: Total enzyme, total substrate pools
• System parameters: Reactor volume, flow rates

Assignment Equations

Create derived quantities that depend on other model components by either using the add_assignment method or the add_assignments method:

# Method 1: Individual assignments
model.add_assignment("E_total", "E + ES + EI")
model.add_assignment("velocity", "(k_cat * ES)")
model.add_assignment("K_eq", "k_forward / k_reverse")

# Method 2: Multiple assignments at once
model.add_assignments(
    E_total="E + ES + EI",
    velocity="(k_cat * ES)",
    K_eq="k_forward / k_reverse"
)

Assignments enable:

• Conservation relationships: Mass balance and stoichiometric constraints
• Derived quantities: Calculated properties for analysis and visualization
• Complex expressions: Multi-parameter relationships and equilibrium definitions

Model Reparameterization

Transform models by substituting symbols with expressions or values:

# Fix specific parameter values
simplified_model = model.reparametrize(k_cat=7.0, K_m=100.0)

# Replace with more complex expressions
regulated_model = model.reparametrize(
    k_cat="k_cat_max * regulator / (K_reg + regulator)"
)

# Eliminate terms by setting to zero
reduced_model = model.reparametrize(inhibition_constant=None)

Reparameterization supports:

• Model simplification: Fixing well-known parameters
• Sensitivity analysis: Systematic parameter variation
• Model comparison: Testing different mechanistic hypotheses

Simulation and Analysis

Basic Simulation Workflow

Run model simulations using the Dataset interface:

from catalax.model.simconfig import SimulationConfig

# Create simulation configuration
config = SimulationConfig(
    t0=0,        # Start time
    t1=100,      # End time  
    nsteps=200,  # Number of time points
    dt0=0.1      # Initial time step
)

# Create dataset with initial conditions
dataset = ctx.Dataset.from_model(model)
dataset.add_initial(S=100.0, E=1.0, P=0.0, ES=0.0)

# Run simulation
results = model.simulate(dataset, config)

# Visualize results
results.plot(show=True)

Model as Predictor

The Model class implements the Predictor interface, enabling use in analysis workflows:

# Generate predictions for new conditions
predictions = model.predict(dataset, config)

# Calculate prediction metrics
metrics = dataset.metrics(model)
print(f"Model R²: {metrics.r2:.3f}")
print(f"Model RMSE: {metrics.rmse:.6f}")

Rate Function Evaluation

Access instantaneous reaction rates for specialized analysis:

import jax.numpy as jnp

# Define time and state points
times = jnp.linspace(0, 100, 50)
states = jnp.array([[100.0, 1.0, 0.0, 0.0]])  # [S, E, P, ES]

# Calculate rates at specific points
rates = model.rates(times, states)
print(f"Rate dimensions: {rates.shape}")  # [n_times, n_states]

Rate evaluation enables:

• Phase plane analysis: Vector field visualization
• Stability analysis: Equilibrium point characterization
• Sensitivity analysis: Gradient-based parameter studies

Model Persistence and Exchange

Saving and Loading Models

Preserve models for reproducible research and collaboration:

# Save model to JSON format
model.save("./models/", "enzyme_kinetics_v1")

# Load model from file
loaded_model = ctx.Model.load("./models/enzyme_kinetics_v1.json")

# Verify model integrity
assert loaded_model.name == model.name
assert len(loaded_model.parameters) == len(model.parameters)

EnzymeML Integration

Exchange models with the broader biochemical modeling community by using the EnzymeML format. You can either create a new EnzymeML document from the model or update an existing EnzymeML document with the model’s parameters:

Method 1: Create a new EnzymeML document from the model

import pyenzyme as pe

# Create EnzymeML document from model
enzmldoc = model.to_enzymeml()

# Load model from existing EnzymeML document
enzymeml_model = ctx.Model.from_enzymeml(enzmldoc)

Method 2: Update an existing EnzymeML document with the model’s parameters

# Assume we have an existing EnzymeML document
enzmldoc = pe.EnzymeMLDocument(name="Enzyme Study")

# Update the existing EnzymeML document with the model's contents
model.update_enzymeml_parameters(enzmldoc)

The EnzymeML integration supports:

• Standards compliance: FAIR data principles and community standards
• Tool interoperability: Exchange with other modeling platforms
• Metadata preservation: Experimental context and provenance information

Best Practices and Guidelines

Follow systematic approaches for robust model development:

1. Start simple: Begin with minimal mechanisms and add complexity incrementally
2. Validate incrementally: Test each addition with simulations and basic analysis
3. Document assumptions: Use clear names and comments to explain modeling choices
4. Version control: Save model states at key development milestones
5. Test reproducibility: Verify that saved/loaded models produce identical results