Model-Based Design of Experiments

Model-Based Design of Experiments#

discopt.doe provides optimal design of experiments for models with unknown parameters, using exact JAX autodiff for Fisher Information Matrix computation. It builds on the parameter-estimation layer in discopt.estimate (part of the base discopt package): both share the same Experiment interface, so a model defined once serves estimation and design alike.

The Experiment Interface (recap)#

You subclass discopt.estimate.Experiment and implement create_model(), which returns an ExperimentModel with labeled components:

Component	Type	Purpose
`unknown_parameters`	`dict[str, Variable]`	Parameters to estimate (optimized as variables)
`design_inputs`	`dict[str, Variable]`	Experimental conditions the experimenter controls
`responses`	`dict[str, Expression]`	Model predictions at measurement points
`measurement_error`	`dict[str, float]`	Standard deviation \(\sigma_i\) for each response

from discopt.estimate import Experiment, ExperimentModel
import discopt.modeling as dm

class MyExperiment(Experiment):
    def create_model(self, **kwargs):
        m = dm.Model("my_exp")
        k = m.continuous("k", lb=0, ub=10)       # unknown parameter
        T = m.continuous("T", lb=300, ub=400)     # design input
        y_pred = k * dm.exp(-1000 / T)            # model prediction
        return ExperimentModel(
            model=m,
            unknown_parameters={"k": k},
            design_inputs={"T": T},
            responses={"y": y_pred},
            measurement_error={"y": 0.05},
        )

See the base discopt documentation for the estimation workflow (estimate_parameters, confidence intervals, covariance).

Fisher Information Matrix#

The FIM quantifies how much information an experiment provides about unknown parameters [Franceschini and Macchietto, 2008]:

\[\text{FIM} = J^T \Sigma^{-1} J\]

where \(J_{ij} = \partial y_i / \partial \theta_j\) is the sensitivity Jacobian. discopt computes \(J\) via exact JAX autodiff — no finite differences, no step-size tuning, no extra model solves [Wang and Dowling, 2022].

from discopt.doe import compute_fim

fim_result = compute_fim(experiment, {"k": 2.0}, {"T": 350.0})
print(fim_result.d_optimal)     # log det(FIM) — D-optimality
print(fim_result.a_optimal)     # trace(FIM^{-1}) — A-optimality
print(fim_result.e_optimal)     # min eigenvalue — E-optimality

Design Criteria [Atkinson et al., 2007]#

Criterion	Formula	Interpretation
D-optimal	\(\max \log \det(\text{FIM})\)	Minimize volume of confidence ellipsoid
A-optimal	\(\min \text{tr}(\text{FIM}^{-1})\)	Minimize average parameter variance
E-optimal	\(\max \lambda_{\min}(\text{FIM})\)	Minimize worst-case variance
ME-optimal	\(\min \kappa(\text{FIM})\)	Balance information across parameters

Optimal Experimental Design#

optimal_experiment() finds design conditions that maximize the chosen information criterion:

from discopt.doe import optimal_experiment, DesignCriterion

design = optimal_experiment(
    experiment,
    param_values={"k": 2.0},
    design_bounds={"T": (300, 400)},
    criterion=DesignCriterion.D_OPTIMAL,
)
print(design.design)                    # {"T": 387.2}
print(design.predicted_standard_errors) # predicted SE if experiment is run

Design Space Exploration#

Visualize how FIM metrics vary across the design space:

from discopt.doe import explore_design_space
import numpy as np

result = explore_design_space(
    experiment,
    param_values={"k": 2.0},
    design_ranges={"T": np.linspace(300, 400, 20)},
)
result.plot_sensitivity(criterion="log_det_fim")

Sequential DoE#

The most powerful workflow alternates estimation and design:

Estimate parameters from all collected data
Compute FIM at current estimates
Design the next experiment to maximize information gain
Run the experiment and collect new data
Repeat — confidence intervals shrink each round

from discopt.doe import sequential_doe

history = sequential_doe(
    experiment=exp,
    initial_data=data,
    initial_guess={"k": 1.0},
    design_bounds={"T": (300, 400)},
    n_rounds=5,
    run_experiment=my_lab_runner,  # callable that runs real experiments
)
for r in history:
    print(f"Round {r.round}: k={r.estimation.parameters['k']:.3f}")

Identifiability Analysis#

Before running experiments, check that parameters are structurally identifiable from the proposed measurements:

from discopt.doe import check_identifiability

result = check_identifiability(experiment, {"k": 2.0, "A": 5.0})
print(result.is_identifiable)          # True/False
print(result.fim_rank)                 # rank of FIM
print(result.problematic_parameters)   # unidentifiable params

API Reference#

Key entry points:

discopt.doe

compute_fim() — Fisher Information Matrix via JAX autodiff
optimal_experiment() — find best experimental conditions
explore_design_space() — grid evaluation with plotting
sequential_doe() — full estimate-design loop
check_identifiability() — parameter identifiability check
DesignCriterion — D/A/E/ME optimality constants
FIMResult — FIM with all optimality metrics

discopt.estimate (base package)

Experiment — base class for experiment definitions
ExperimentModel — annotated model with metadata
estimate_parameters() — run parameter estimation

Comparison with Pyomo.DoE#

Feature	Pyomo.DoE [Wang and Dowling, 2022]	discopt.doe
Sensitivity Jacobian	Finite differences (2N extra solves)	Exact JAX autodiff
FIM gradient w.r.t. design	Cholesky variables in NLP	JAX autodiff through FIM
GPU acceleration	No	Yes (JAX backend)
Parameter estimation	Separate package (parmest)	Integrated (`discopt.estimate`)
Mixed-integer designs	Not supported	Native MINLP
Implicit differentiation	Not available	Envelope theorem