Introduction

odesign is an optimal design of experiments library written in pure rust. It allows to find optimal designs for (quasi-) linear models considering arbitrary optimalities.

This book serves as a high level introduction and theoretical background.

Please find more resources here:

docs on docs.rs/odesign
project on sourcehut sr.ht/~maunke/odesign
read-only mirror on github maunke/odesign

Use Cases

Fast calculation of optimal designs of arbitrary linear models with custom design bounds and optimalities.
Research in area of optimal designs; e.g. I am working on a new optimal design feature selection algorithm, a mixture of SFFS, D-, C- and Measurements-Costs-Optimality, allowing to perform model feature selection and measurements alternating.

Core Features

The library consists of three main features:

Basic example

In short, this is a basic example of an optimal design of the simple polynomial 1 + x within design bounds [-1, +1] and 101 equally distributed grid points as an init design.

use nalgebra::{SVector, Vector1};
use num_dual::DualNum;
use odesign::{
    DOptimality, Feature, FeatureFunction, FeatureSet, LinearModel, OptimalDesign, Result,
};
use std::sync::Arc;

#[derive(Feature)]
#[dimension = 1]
struct Monomial {
    i: i32,
}

impl FeatureFunction<1> for Monomial {
    fn f<D: DualNum<f64>>(&self, x: &SVector<D, 1>) -> D {
        x[0].powi(self.i)
    }
}

// f(x): 1 + x
fn main() -> Result<()> {
    let mut fs = FeatureSet::new();
    let c: Arc<_> = Monomial { i: 0 }.into();
    fs.push(c);
    let c: Arc<_> = Monomial { i: 1 }.into();
    fs.push(c);

    let lm = LinearModel::new(fs.features);

    let optimality: Arc<_> = DOptimality::new(lm.into()).into();
    let lower = Vector1::new(-1.0);
    let upper = Vector1::new(1.0);
    let q = Vector1::new(101);

    let mut od = OptimalDesign::new()
        .with_optimality(optimality)
        .with_bound_args(lower, upper)?
        .with_init_design_grid_args(lower, upper, q)?;
    od.solve();

    println!("{od}");

    Ok(())
}
// Output
// ---------- Design ----------
// Weight  Support Vector
// 0.5000  [ -1.0000 ]
// 0.5000  [ +1.0000 ]
// -------- Statistics --------
// Optimality measure: 1.000000
// No. support vectors: 2
// Iterations: 1
// ----------------------------

Feature derive

In order to define a differentiable linear model with its features you just need to define a set of instantiated structs with the feature derive and implement their feature function. With help of this design arbitrary linear models can be defined. See here the simple implementation of the following linear model

$$ f(x) = \beta_0 + \beta_1 \cdot x + \beta_2 \cdot \frac{1}{x}$$

with its coefficient $\beta$ we want to estimate in the context of linear regression.

#[derive(Feature)]
#[dimension = 1]
struct Monomial {
    i: i32,
}

impl FeatureFunction<1> for Monomial {
    fn f<D: DualNum<f64>>(&self, x: &SVector<D, 1>) -> D {
        x[0].powi(self.i)
    }
}

/// f(x): 1 + x + 1 / x
fn main() -> Result<()> {
    let mut fs = FeatureSet::new();
    for i in -1..2 {
        let c: Arc<_> = Monomial { i }.into();
        fs.push(c);
    }

    let _lm = LinearModel::new(fs.features);

    Ok(())
}

Optimalities

Matrix Means

The best known class of information functions was introduced by Kiefer, 1974 (doi: 10.1214/aos/1176342810). Each information function is described by a hyperparameter δ ∈ (-∞, 1], which can be categorized as follows into the class of matrix means ∆δ : PD (n) → R+ (where PD(n) stands for a positive definite symmetric n x n matrix) as follows

$$ \Delta_{\delta}\left(\mathcal{M}(\xi)\right) = \left[\frac{1}{n}\text{tr}{\mathcal{M}(\xi)^\delta}\right]^{\frac{1}{\delta}} \quad \text{,} $$

with the fisher information matrix M and the design ξ.

Right now the popular D-Optimalty with δ = 0 is part of odesign. Additionally the C-Optimality with

$$ \Delta_c(\xi) := c^T \mathcal{M}^{-1}(\xi) c $$

is implemented.

Custom optimality

Beside the predefined optimalities, you can create your custom ones and/or create a weighted sum of different optimalities, e.g. CD-Optimality.

Lets consider a simple example where we define a to be maximized costs-efficiency-optimality that is concave and sums up the negative weighted sum of the norm of each support vector

$$ \Delta_{\text{costs}} := \exp\big(-\sum_{i \in [1, .., n]} w_i \cdot || x_i || \big) $$

with n support vectors x and theirs weights w (see the custom-optimality example here) or run the following command:

cargo run --example custom-optimality --release

Since the solver minimizes the sum of negative log of the desired optimalities, we formulate the custom matrix means derivatives as follows:

Value:

$$ -\log \Delta_{\text{costs}} = \sum_{i \in [1, .., n]} w_i \cdot || x_i || $$

Gradient:

$$ \frac{\partial -\log \Delta_{\text{costs}}}{\partial w_i} = || x_i || $$

Hessian:

$$ \frac{\partial^2 -\log \Delta_{\text{costs}}}{\partial w_i \partial w_j} = 0 $$

The same principle applies to the dispersion function, where we will derive to the support vector.

Optimal Design Solver

This library implements the "Adaptive grid Semidefinite Programming for finding optimal designs" method (doi: 10.1007/s11222-017-9741-y), please find the open accessible paper here.

odesign Book